Packages
if (!require("tidyverse")) install.packages("tidyverse")
library("tidyverse") # for tidy/dplyr workflow of calculations
if (!require("lmerTest")) install.packages("lmerTest")
library("lmerTest") # LMMs
if (!require("ggplot2")) install.packages("ggplot2")
library("ggplot2") # plots
if (!require("RColorBrewer")) install.packages("RColorBrewer")
library("RColorBrewer") # pretty colors
if (!require("ggpubr")) install.packages("ggpubr")
library("ggpubr") # for multi-ggplot figs
if (!require("emmeans")) install.packages("emmeans")
library("emmeans") # marginal means & confintsSet Colors & Themes
savs_ggplot_theme <- theme_classic() +
theme(text = element_text(color = "black",
family = "sans",
size = 12),
axis.text = element_text(color = "black",
family = "sans",
size = 8),
legend.text = element_text(color = "black",
family = "sans",
size = 8),
legend.text.align = 0,
legend.position = "bottom"
)
# custom colors
Male <- brewer.pal(11, "RdYlBu")[c(10)]
Fem <- brewer.pal(11, "RdYlBu")[c(3)]
Fem_grav <- brewer.pal(11, "RdYlBu")[c(2)]
Fem_nongrav <- brewer.pal(11, "RdYlBu")[c(4)]
my_colors <- c(Fem, Male)
my_shapes <- c(19,15)
my_shapes_open <- c(21,22)
# corr colors
n3 <- brewer.pal(11, "RdYlBu")[c(9, 11, 3)]
n2 <- brewer.pal(9, "YlGnBu")[c(6, 8)]
n1 <- brewer.pal(9, "YlGnBu")[c(6)]
# labels for plots for repeat measures lizards only
my_labels <- c("April\nF\n(n=11)", "April\nM\n(n=22)",
"May\nF\n(n=10)", "May\nM\n(n=17)",
"July\nF\n(n=2)", "July\nM\n(n=12)")
# labels for plots for full dataset of all lizards
my_labels_full <- c("April\nF\n(n=26)", "April\nM\n(n=41)",
"May\nF\n(n=15)", "May\nM\n(n=22)",
"July\nF\n(n=3)", "July\nM\n(n=12)")
my_labels_full_F <- c("April\n(n=26)",
"May\n(n=15)",
"July\n(n=3)")
my_labels_full_M <- c("April\n(n=41)",
"May\n(n=22)",
"July\n(n=12)")Background & Introduction
This data was collected on Blunt-nosed Leopard Lizards (Gambelia sila) during their active season in 2021. In this document, we will investigate the nuances and trends in the hydric physiology of these desert lizards.
We will do some additional data wrangling (see other ‘wrangling_’ files in this repository for more details on data preparation), compute some initial summary statistics, and investigate the following questions:
- Are the different hydric physiology metrics we measured related to each other?
- Did our attempt at supplemental hydration actually do anything?
- How did hydric physiology change over time throughout the active season (April-July) of these lizards? This section composes the biggest chunk of this file.
- Does hydric physiology relate to reproductive ability or effort in any way?
Load in Data
Main Dataframe
dat <- read_rds("./data/G_sila_clean_full_dat.RDS") %>%
mutate(month_sex = factor(paste(month, sex_M_F),
levels = c("April Female", "April Male",
"May Female", "May Male",
"July Female", "July Male")),
date_chunk_pretty = factor(case_when(month == "April" ~ "Apr 23-25",
month == "May" ~ "May 7-9",
month == "July" ~ "Jul 14"),
levels = c("Apr 23-25", "May 7-9", "Jul 14")),
gravid_Y_N = factor(gravid_Y_N, levels = c("Not Gravid", "Gravid")),
percent_water_mass_change = (tmt_mass_change_g/mass_g)*100)
summary(dat)## capture_date_time capture_date radio_collar_serial
## Min. :2021-04-23 03:39:00.00 Min. :2021-04-23 Length:119
## 1st Qu.:2021-04-23 08:39:15.00 1st Qu.:2021-04-24 Class :character
## Median :2021-04-24 05:22:00.00 Median :2021-04-24 Mode :character
## Mean :2021-04-28 07:08:18.98 Mean :2021-05-08
## 3rd Qu.:2021-05-07 06:16:45.00 3rd Qu.:2021-05-08
## Max. :2021-05-08 07:19:00.00 Max. :2021-07-14
## NA's :21
## PIT_tag_ID individual_ID sex_M_F mass_g SVL_mm
## Length:119 F-12 : 3 Female:44 Min. :20.00 Min. : 77.0
## Class :character M-09 : 3 Male :75 1st Qu.:31.00 1st Qu.: 93.5
## Mode :character M-10 : 3 Median :37.00 Median : 99.0
## M-11 : 3 Mean :36.94 Mean : 99.0
## M-19 : 3 3rd Qu.:43.00 3rd Qu.:106.5
## M-20 : 3 Max. :58.80 Max. :122.0
## (Other):101
## tmt tmt_mass_change_g capture_to_msmt hematocrit_percent
## Water Tmt:38 Min. :-3.0000 Min. : 10.00 Min. :10.00
## Sham Tmt :38 1st Qu.:-0.1500 1st Qu.: 49.25 1st Qu.:30.00
## No Tmt :43 Median : 0.0000 Median : 117.50 Median :34.00
## Mean :-0.1327 Mean : 276.18 Mean :33.78
## 3rd Qu.: 0.0000 3rd Qu.: 228.75 3rd Qu.:38.00
## Max. : 3.0000 Max. :1665.00 Max. :58.00
## NA's :15 NA's :21 NA's :3
## cloacal_temp_C CEWL_g_m2h msmt_temp_C msmt_RH_percent
## Min. :22.00 Min. : 0.650 Min. :19.00 Min. :11.84
## 1st Qu.:31.00 1st Qu.: 8.486 1st Qu.:28.57 1st Qu.:14.64
## Median :33.00 Median :10.381 Median :30.38 Median :16.12
## Mean :32.22 Mean :10.272 Mean :29.43 Mean :20.90
## 3rd Qu.:34.00 3rd Qu.:12.881 3rd Qu.:31.52 3rd Qu.:23.78
## Max. :36.00 Max. :21.673 Max. :33.55 Max. :40.22
## NA's :1 NA's :1 NA's :1
## e_s_kPa msmt_VPD_kPa osmolality_mmol_kg_mean month week
## Min. :2.196 Min. :1.371 Min. :309.3 April:67 1 :67
## 1st Qu.:3.906 1st Qu.:2.946 1st Qu.:348.2 May :37 3 :37
## Median :4.336 Median :3.635 Median :366.5 July :15 12:15
## Mean :4.160 Mean :3.332 Mean :367.8
## 3rd Qu.:4.626 3rd Qu.:3.960 3rd Qu.:379.9
## Max. :5.188 Max. :4.319 Max. :437.3
## NA's :1 NA's :1 NA's :3
## week_num ultrasound_date gravid_Y_N number_eggs
## Min. : 1.000 Min. :2021-04-24 Not Gravid:12 Min. :3.000
## 1st Qu.: 1.000 1st Qu.:2021-04-24 Gravid :16 1st Qu.:3.750
## Median : 1.000 Median :2021-04-24 NA's :91 Median :4.000
## Mean : 3.008 Mean :2021-04-29 Mean :3.812
## 3rd Qu.: 3.000 3rd Qu.:2021-05-08 3rd Qu.:4.000
## Max. :12.000 Max. :2021-05-08 Max. :5.000
## NA's :91 NA's :103
## largest_egg_diameter_cm largest_egg_circumference_cm stage
## Min. :0.5700 Min. :1.710 small round: 1
## 1st Qu.:0.8425 1st Qu.:2.385 large round: 9
## Median :0.9950 Median :2.885 soft : 2
## Mean :1.1131 Mean :3.148 soft oblong: 1
## 3rd Qu.:1.4200 3rd Qu.:3.985 firm oblong: 2
## Max. :1.8000 Max. :4.740 hard oblong: 1
## NA's :103 NA's :103 NA's :103
## dev_point SMI month_sex date_chunk_pretty
## Min. :1.00 Min. :22.35 April Female:26 Apr 23-25:67
## 1st Qu.:2.00 1st Qu.:31.58 April Male :41 May 7-9 :37
## Median :2.00 Median :34.29 May Female :15 Jul 14 :15
## Mean :2.75 Mean :34.57 May Male :22
## 3rd Qu.:3.25 3rd Qu.:37.33 July Female : 3
## Max. :5.00 Max. :45.63 July Male :12
## NA's :103
## percent_water_mass_change
## Min. :-6.8337
## 1st Qu.:-0.4487
## Median : 0.0000
## Mean :-0.3281
## 3rd Qu.: 0.0000
## Max. :10.0000
## NA's :15## [1] "245-761" "252-886" "245-762" "252-887" "245-753" "245-752"
## [7] "245-754" "229-079" "245-744" "229-069" "252-885" "245-759"
## [13] "245-751" "245-757" "229-059" "229-063" "245-747" "229-070"
## [19] "229-065" "245-741" "245-748" "252-882" "252-881" "229-078"
## [25] "245-750" "229-080" "245-760" "245-746" "245-745" "245-742"
## [31] "252-884" "252-883" "245-756" "229-073" "229-066" "229-072"
## [37] "229-060" "229-074" "229-077" "245-762a" "229-070a" NAdat_F <- dat %>%
dplyr::filter(sex_M_F == "Female")
dat_M <- dat %>%
dplyr::filter(sex_M_F == "Male")
# summary statistics
dat %>%
summarise(mean(CEWL_g_m2h, na.rm = T), sd(CEWL_g_m2h, na.rm = T),
mean(osmolality_mmol_kg_mean, na.rm = T),
sd(osmolality_mmol_kg_mean, na.rm = T),
mean(hematocrit_percent, na.rm = T), sd(hematocrit_percent, na.rm = T),
mean(mass_g), sd(mass_g),
mean(SVL_mm), sd(SVL_mm),
mean(SMI), sd(SMI))## mean(CEWL_g_m2h, na.rm = T) sd(CEWL_g_m2h, na.rm = T)
## 1 10.27202 3.924563
## mean(osmolality_mmol_kg_mean, na.rm = T)
## 1 367.8477
## sd(osmolality_mmol_kg_mean, na.rm = T) mean(hematocrit_percent, na.rm = T)
## 1 26.73921 33.77586
## sd(hematocrit_percent, na.rm = T) mean(mass_g) sd(mass_g) mean(SVL_mm)
## 1 7.363938 36.94034 8.972626 99
## sd(SVL_mm) mean(SMI) sd(SMI)
## 1 9.401587 34.56663 4.475429Variable Explanation
This is an explanation of each variable in our dataframe:
- “capture_date_time” = date and time the lizard was captured
- “capture_date” = date of capture
- “radio_collar_serial” = serial number of the radio-collar that lizard was fitted with, if applicable; only some lizards were fitted with collars
- “PIT_tag_ID” = serial number of the PIT tag injected into the lizard at initial capture; all lizards got PIT-tagged
- “individual_ID” = ID number assigned to lizards
- “sex_M_F” = lizard sex
- “mass_g” = lizard body mass
- “SVL_mm” = lizard snout-vent length
- “tmt” = treatment given to the lizard: water supplementation, sham, or none (for the lizards that were not radio-collared)
- “tmt_mass_change_g” = the change in mass for lizards from before to after treatment; only for lizards given water or sham treatments
- “capture_to_msmt” = time between a lizard being captured and being measured
- “hematocrit_percent” = percent of blood sample that was red blood cells
- “cloacal_temp_C” = internal body temperature of the lizard at the time CEWL was taken, taken with a thermocouple in the cloaca
- “CEWL_g_m2h” = cutaneous evaporative water loss; see wrangling_CEWL for data cleaning
- “msmt_temp_C” = ambient temperature at the time of CEWL
measurement
- “msmt_RH_percent” = ambient humidity at the time of CEWL
measurement
- “e_s_kPa” = saturation vapor pressure for that ambient temperature
- “msmt_VPD_kPa” = ambient vapor pressure deficit at the time of CEWL
measurement
- “osmolality_mmol_kg_mean” = plasma osmolality; see wrangling_osml for data cleaning
- “month” = measurement period as categorical month (April, May, June)
- “week” = measurement period as categorical week throughout the summer (1, 3, 12)
- “week_num” = measurement period as numerical week
- “ultrasound_date” = date that a female was ultrasounded; may differ from capture and measurement dates
- “gravid_Y_N” = whether or not a female was measured with eggs, if
measured (not measured for the non-radio-collared ones)
- “number_eggs” = if female is gravid, the number of eggs visible with the ultrasound
- “largest_egg_diameter_cm” = if female is gravid, the diameter of the largest egg in the clutch
- “largest_egg_circumference_cm” = if female is gravid, the circumference of the largest egg in the clutch
- “stage” = if female is gravid, approximate stage of development of the eggs
- “dev_point” = numerical version of development stage
- “SMI” = body condition of lizards as scaled mass index (calculated in wrangling_general)
- “month_sex” = variable that I created to make categories for each sex each month
- “date_chunk_pretty” = nicely formatted measurement period labels
Repeat Subset
I want to make some sub-dataframes for when I want to do stats for only a subset of individuals.
First, get a dataframe of only individuals that had more than one measurement taken, and note how many individuals had each n measurements:
# how many lizards were measured 1, 2, or 3 times?
note_repeats <- dat %>%
group_by(individual_ID) %>%
summarise(n = n()) %>%
group_by(n) %>%
summarise(n_lizards = n())
note_repeats ## # A tibble: 3 × 2
## n n_lizards
## <int> <int>
## 1 1 45
## 2 2 28
## 3 3 6## [1] 79# get subsetted data
dat_repeats <- dat %>%
group_by(individual_ID) %>%
mutate(n = n()) %>%
dplyr::filter(n>1)
# check these match 'note_repeats'
dat_repeats %>%
group_by(individual_ID) %>%
summarise(n = n()) %>%
group_by(n) %>%
summarise(n_lizards = n())## # A tibble: 2 × 2
## n n_lizards
## <int> <int>
## 1 2 28
## 2 3 6Female Subset
To look at reproductive effort ~ water balance, I will make another subset for only the females we collected gravidity/egg data on:
dat_repro_females <- dat %>%
dplyr::filter(sex_M_F == "Female") %>%
dplyr::filter(complete.cases(gravid_Y_N)) %>%
group_by(individual_ID) %>%
mutate(n = n(),
n = factor(n))
summary(dat_repro_females)## capture_date_time capture_date radio_collar_serial
## Min. :2021-04-23 04:19:00.0 Min. :2021-04-23 Length:28
## 1st Qu.:2021-04-23 08:30:00.0 1st Qu.:2021-04-23 Class :character
## Median :2021-04-24 04:11:00.0 Median :2021-04-24 Mode :character
## Mean :2021-04-28 09:34:38.4 Mean :2021-04-29
## 3rd Qu.:2021-05-08 01:54:00.0 3rd Qu.:2021-05-08
## Max. :2021-05-08 07:19:00.0 Max. :2021-05-08
## NA's :3
## PIT_tag_ID individual_ID sex_M_F mass_g SVL_mm
## Length:28 F-01 : 2 Female:28 Min. :26.00 Min. : 85.0
## Class :character F-03 : 2 Male : 0 1st Qu.:33.52 1st Qu.: 95.0
## Mode :character F-05 : 2 Median :38.00 Median :100.0
## F-06 : 2 Mean :39.79 Mean :101.1
## F-11 : 2 3rd Qu.:44.50 3rd Qu.:108.0
## F-12 : 2 Max. :58.80 Max. :114.0
## (Other):16
## tmt tmt_mass_change_g capture_to_msmt hematocrit_percent
## Water Tmt:15 Min. :-3.0000 Min. : 25.0 Min. :18.00
## Sham Tmt :13 1st Qu.:-1.0000 1st Qu.: 49.0 1st Qu.:26.75
## No Tmt : 0 Median :-0.1000 Median : 118.0 Median :32.00
## Mean :-0.4519 Mean : 238.8 Mean :32.14
## 3rd Qu.: 0.0000 3rd Qu.: 163.0 3rd Qu.:37.00
## Max. : 1.9000 Max. :1665.0 Max. :44.00
## NA's :1 NA's :3
## cloacal_temp_C CEWL_g_m2h msmt_temp_C msmt_RH_percent
## Min. :28.00 Min. : 1.493 Min. :23.13 Min. :12.20
## 1st Qu.:33.00 1st Qu.: 9.171 1st Qu.:30.36 1st Qu.:14.38
## Median :34.00 Median :10.239 Median :31.40 Median :15.55
## Mean :33.46 Mean :11.423 Mean :30.70 Mean :16.98
## 3rd Qu.:35.00 3rd Qu.:13.613 3rd Qu.:31.93 3rd Qu.:19.05
## Max. :36.00 Max. :21.673 Max. :33.55 Max. :28.17
##
## e_s_kPa msmt_VPD_kPa osmolality_mmol_kg_mean month week
## Min. :2.830 Min. :2.033 Min. :309.3 April:17 1 :17
## 1st Qu.:4.331 1st Qu.:3.600 1st Qu.:337.8 May :11 3 :11
## Median :4.595 Median :3.863 Median :346.6 July : 0 12: 0
## Mean :4.444 Mean :3.705 Mean :352.6
## 3rd Qu.:4.736 3rd Qu.:4.050 3rd Qu.:369.2
## Max. :5.188 Max. :4.319 Max. :402.0
##
## week_num ultrasound_date gravid_Y_N number_eggs
## Min. :1.000 Min. :2021-04-24 Not Gravid:12 Min. :3.000
## 1st Qu.:1.000 1st Qu.:2021-04-24 Gravid :16 1st Qu.:3.750
## Median :1.000 Median :2021-04-24 Median :4.000
## Mean :1.786 Mean :2021-04-29 Mean :3.812
## 3rd Qu.:3.000 3rd Qu.:2021-05-08 3rd Qu.:4.000
## Max. :3.000 Max. :2021-05-08 Max. :5.000
## NA's :12
## largest_egg_diameter_cm largest_egg_circumference_cm stage
## Min. :0.5700 Min. :1.710 small round: 1
## 1st Qu.:0.8425 1st Qu.:2.385 large round: 9
## Median :0.9950 Median :2.885 soft : 2
## Mean :1.1131 Mean :3.148 soft oblong: 1
## 3rd Qu.:1.4200 3rd Qu.:3.985 firm oblong: 2
## Max. :1.8000 Max. :4.740 hard oblong: 1
## NA's :12 NA's :12 NA's :12
## dev_point SMI month_sex date_chunk_pretty
## Min. :1.00 Min. :27.04 April Female:17 Apr 23-25:17
## 1st Qu.:2.00 1st Qu.:32.45 April Male : 0 May 7-9 :11
## Median :2.00 Median :34.89 May Female :11 Jul 14 : 0
## Mean :2.75 Mean :35.82 May Male : 0
## 3rd Qu.:3.25 3rd Qu.:40.19 July Female : 0
## Max. :5.00 Max. :45.63 July Male : 0
## NA's :12
## percent_water_mass_change n
## Min. :-6.834 2:18
## 1st Qu.:-3.046 1:10
## Median :-0.250
## Mean :-1.171
## 3rd Qu.: 0.000
## Max. : 7.037
## NA's :1It’s good to note that I measured reproductive effort at two time points for 9 females, and only at one time point for 10 females.
Check whether the single-measurement females are all from the first weekend or not:
##
## 2 1
## 1 9 8
## 3 9 2
## 12 0 0## capture_date_time capture_date radio_collar_serial
## Min. :2021-04-23 04:19:00 Min. :2021-04-23 Length:10
## 1st Qu.:2021-04-24 02:57:00 1st Qu.:2021-04-24 Class :character
## Median :2021-04-24 04:10:00 Median :2021-04-24 Mode :character
## Mean :2021-04-25 13:16:20 Mean :2021-04-26
## 3rd Qu.:2021-04-24 05:56:00 3rd Qu.:2021-04-24
## Max. :2021-05-08 07:19:00 Max. :2021-05-08
## NA's :1
## PIT_tag_ID individual_ID sex_M_F mass_g SVL_mm
## Length:10 F-02 :1 Female:10 Min. :26.00 Min. : 85.00
## Class :character F-04 :1 Male : 0 1st Qu.:32.00 1st Qu.: 96.25
## Mode :character F-07 :1 Median :37.50 Median :103.50
## F-08 :1 Mean :36.30 Mean :101.80
## F-08a :1 3rd Qu.:41.75 3rd Qu.:110.25
## F-09 :1 Max. :47.00 Max. :114.00
## (Other):4
## tmt tmt_mass_change_g capture_to_msmt hematocrit_percent
## Water Tmt:7 Min. :-1.00 Min. : 35.0 Min. :18.00
## Sham Tmt :3 1st Qu.:-1.00 1st Qu.: 49.0 1st Qu.:25.00
## No Tmt :0 Median :-0.50 Median : 91.0 Median :26.00
## Mean :-0.21 Mean : 261.7 Mean :29.20
## 3rd Qu.: 0.00 3rd Qu.: 132.0 3rd Qu.:35.75
## Max. : 1.90 Max. :1665.0 Max. :41.00
## NA's :1
## cloacal_temp_C CEWL_g_m2h msmt_temp_C msmt_RH_percent
## Min. :28.00 Min. : 8.457 Min. :25.03 Min. :14.12
## 1st Qu.:32.25 1st Qu.: 9.554 1st Qu.:29.50 1st Qu.:15.05
## Median :33.50 Median :10.159 Median :30.96 Median :16.91
## Mean :33.00 Mean :11.986 Mean :30.54 Mean :18.24
## 3rd Qu.:34.00 3rd Qu.:14.551 3rd Qu.:32.12 3rd Qu.:20.05
## Max. :35.00 Max. :18.743 Max. :33.55 Max. :26.03
##
## e_s_kPa msmt_VPD_kPa osmolality_mmol_kg_mean month week
## Min. :3.172 Min. :2.346 Min. :309.3 April:8 1 :8
## 1st Qu.:4.124 1st Qu.:3.341 1st Qu.:330.6 May :2 3 :2
## Median :4.482 Median :3.776 Median :357.3 July :0 12:0
## Mean :4.410 Mean :3.623 Mean :355.6
## 3rd Qu.:4.786 3rd Qu.:4.014 3rd Qu.:378.2
## Max. :5.188 Max. :4.319 Max. :402.0
##
## week_num ultrasound_date gravid_Y_N number_eggs
## Min. :1.0 Min. :2021-04-24 Not Gravid:6 Min. :3.00
## 1st Qu.:1.0 1st Qu.:2021-04-24 Gravid :4 1st Qu.:3.75
## Median :1.0 Median :2021-04-24 Median :4.00
## Mean :1.4 Mean :2021-04-26 Mean :3.75
## 3rd Qu.:1.0 3rd Qu.:2021-04-24 3rd Qu.:4.00
## Max. :3.0 Max. :2021-05-08 Max. :4.00
## NA's :6
## largest_egg_diameter_cm largest_egg_circumference_cm stage
## Min. :0.570 Min. :1.710 small round:0
## 1st Qu.:0.795 1st Qu.:2.265 large round:3
## Median :0.935 Median :2.595 soft :0
## Mean :0.935 Mean :2.717 soft oblong:0
## 3rd Qu.:1.075 3rd Qu.:3.047 firm oblong:1
## Max. :1.300 Max. :3.970 hard oblong:0
## NA's :6 NA's :6 NA's :6
## dev_point SMI month_sex date_chunk_pretty
## Min. :2.00 Min. :27.04 April Female:8 Apr 23-25:8
## 1st Qu.:2.00 1st Qu.:30.81 April Male :0 May 7-9 :2
## Median :2.00 Median :32.93 May Female :2 Jul 14 :0
## Mean :2.75 Mean :32.21 May Male :0
## 3rd Qu.:2.75 3rd Qu.:33.44 July Female :0
## Max. :5.00 Max. :37.68 July Male :0
## NA's :6
## percent_water_mass_change n
## Min. :-3.8462 2: 0
## 1st Qu.:-2.6849 1:10
## Median :-1.1628
## Mean :-0.5466
## 3rd Qu.: 0.0000
## Max. : 7.0370
## ## capture_date_time capture_date radio_collar_serial
## Min. :2021-05-08 01:40:00.00 Min. :2021-05-08 Length:9
## 1st Qu.:2021-05-08 02:09:00.00 1st Qu.:2021-05-08 Class :character
## Median :2021-05-08 02:58:00.00 Median :2021-05-08 Mode :character
## Mean :2021-05-08 03:06:42.85 Mean :2021-05-08
## 3rd Qu.:2021-05-08 03:42:00.00 3rd Qu.:2021-05-08
## Max. :2021-05-08 05:27:00.00 Max. :2021-05-08
## NA's :2
## PIT_tag_ID individual_ID sex_M_F mass_g SVL_mm
## Length:9 F-01 :1 Female:9 Min. :33.70 Min. : 90.00
## Class :character F-03 :1 Male :0 1st Qu.:35.60 1st Qu.: 95.00
## Mode :character F-05 :1 Median :43.90 Median : 98.00
## F-06 :1 Mean :44.02 Mean : 99.56
## F-11 :1 3rd Qu.:52.70 3rd Qu.:105.00
## F-12 :1 Max. :58.80 Max. :110.00
## (Other):3
## tmt tmt_mass_change_g capture_to_msmt hematocrit_percent
## Water Tmt:4 Min. :-3.0 Min. : 59.0 Min. :31.00
## Sham Tmt :5 1st Qu.:-1.5 1st Qu.:158.5 1st Qu.:32.00
## No Tmt :0 Median :-1.0 Median :244.0 Median :37.00
## Mean :-1.0 Mean :201.7 Mean :37.33
## 3rd Qu.: 0.0 3rd Qu.:259.5 3rd Qu.:43.00
## Max. : 1.0 Max. :273.0 Max. :44.00
## NA's :1 NA's :2
## cloacal_temp_C CEWL_g_m2h msmt_temp_C msmt_RH_percent
## Min. :33.00 Min. : 8.318 Min. :30.40 Min. :13.12
## 1st Qu.:34.00 1st Qu.: 9.318 1st Qu.:30.86 1st Qu.:14.44
## Median :35.00 Median :10.286 Median :31.45 Median :14.70
## Mean :34.56 Mean :10.915 Mean :31.38 Mean :14.78
## 3rd Qu.:35.00 3rd Qu.:12.213 3rd Qu.:31.90 3rd Qu.:15.48
## Max. :36.00 Max. :14.378 Max. :32.04 Max. :15.84
##
## e_s_kPa msmt_VPD_kPa osmolality_mmol_kg_mean month week
## Min. :4.340 Min. :3.685 Min. :337.3 April:0 1 :0
## 1st Qu.:4.456 1st Qu.:3.813 1st Qu.:340.5 May :9 3 :9
## Median :4.608 Median :3.939 Median :361.0 July :0 12:0
## Mean :4.591 Mean :3.913 Mean :357.0
## 3rd Qu.:4.728 3rd Qu.:4.028 3rd Qu.:369.0
## Max. :4.765 Max. :4.112 Max. :379.5
##
## week_num ultrasound_date gravid_Y_N number_eggs
## Min. :3 Min. :2021-05-08 Not Gravid:1 Min. :3.00
## 1st Qu.:3 1st Qu.:2021-05-08 Gravid :8 1st Qu.:3.00
## Median :3 Median :2021-05-08 Median :4.00
## Mean :3 Mean :2021-05-08 Mean :3.75
## 3rd Qu.:3 3rd Qu.:2021-05-08 3rd Qu.:4.00
## Max. :3 Max. :2021-05-08 Max. :5.00
## NA's :1
## largest_egg_diameter_cm largest_egg_circumference_cm stage
## Min. :0.640 Min. :1.820 small round:1
## 1st Qu.:1.135 1st Qu.:3.333 large round:2
## Median :1.440 Median :4.000 soft :2
## Mean :1.335 Mean :3.712 soft oblong:1
## 3rd Qu.:1.560 3rd Qu.:4.228 firm oblong:1
## Max. :1.800 Max. :4.740 hard oblong:1
## NA's :1 NA's :1 NA's :1
## dev_point SMI month_sex date_chunk_pretty
## Min. :1.000 Min. :34.49 April Female:0 Apr 23-25:0
## 1st Qu.:2.000 1st Qu.:36.44 April Male :0 May 7-9 :9
## Median :3.000 Median :41.85 May Female :9 Jul 14 :0
## Mean :3.125 Mean :40.73 May Male :0
## 3rd Qu.:4.250 3rd Qu.:43.87 July Female :0
## Max. :5.000 Max. :45.63 July Male :0
## NA's :1
## percent_water_mass_change n
## Min. :-6.834 2:9
## 1st Qu.:-3.649 1:0
## Median :-2.468
## Mean :-2.203
## 3rd Qu.: 0.000
## Max. : 2.809
## NA's :1##
## 2 1
## Not Gravid 6 6
## Gravid 12 4For the females who were only measured once, 8 of them are from the first weekend and 2 of them are from the third weekend.
Also save a subset of the females that were measured twice:
Does SVL change over time?
I want to know whether SVL changed over time for the individuals who had more than one measurement taken.
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: SVL_mm ~ week_num + (1 | individual_ID)
## Data: dat_repeats
##
## REML criterion at convergence: 495.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.46209 -0.38905 -0.04905 0.51484 2.46481
##
## Random effects:
## Groups Name Variance Std.Dev.
## individual_ID (Intercept) 55.65 7.460
## Residual 19.96 4.467
## Number of obs: 74, groups: individual_ID, 34
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 100.5373 1.4817 42.0626 67.852 <2e-16 ***
## week_num 0.1016 0.1442 42.5756 0.705 0.485
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## week_num -0.359There is no meaningful change in SVL for the individuals we were able to measure, so I could use either mass or body condition to estimate change in body size over time.
Check CEWL Covariates
How much is CEWL impacted by the body temperature, ambient temperature, and ambient VPD at the time of CEWL measurement? I will investigate this using all measurements.
It’s possible that I may find a relationship between ambient temp and CEWL but not body temp and CEWL because the body temp measurements we took were unfortunately only to the nearest degree (vs tenth of a degree).
Plot
First, make a simple plot of each to look for relationships visually.
ggplot(dat) +
aes(x = cloacal_temp_C,
y = CEWL_g_m2h,
color = week
#color = sex_M_F
) +
geom_point() +
geom_smooth(method = "lm") +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'## Warning: Removed 1 rows containing non-finite values (`stat_smooth()`).## Warning: Removed 1 rows containing missing values (`geom_point()`).
ggplot(dat) +
aes(x = msmt_temp_C,
y = CEWL_g_m2h,
color = week
#color = sex_M_F
) +
geom_point() +
geom_smooth(method = "lm") +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'## Warning: Removed 1 rows containing non-finite values (`stat_smooth()`).
## Removed 1 rows containing missing values (`geom_point()`).
ggplot(dat) +
aes(x = msmt_VPD_kPa,
y = CEWL_g_m2h,
color = week
#color = sex_M_F
) +
geom_point() +
geom_smooth(method = "lm") +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'## Warning: Removed 1 rows containing non-finite values (`stat_smooth()`).
## Removed 1 rows containing missing values (`geom_point()`).
Across dates, CEWL had no relationship with body temp, ambient temp, or ambient VPD at the time of capture. However, when separated by date, there were relationships. CEWL had a positive relationship with body temp, ambient temp, and ambient VPD, but only for week 3 measurements. There were no relationships when separated by treatment group.
Models
Now, check linear models of the effect of each of body temp, ambient temp, and ambient VPD with the interaction of date on CEWL.
CEWL_cov_check_clotemp <- lm(data = dat,
CEWL_g_m2h ~ cloacal_temp_C*week)
summary(CEWL_cov_check_clotemp)##
## Call:
## lm(formula = CEWL_g_m2h ~ cloacal_temp_C * week, data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.1993 -2.0150 0.1478 2.0471 10.8241
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 11.45021 4.57266 2.504 0.013718 *
## cloacal_temp_C -0.01821 0.14195 -0.128 0.898154
## week3 -63.35695 17.24290 -3.674 0.000367 ***
## week12 -20.29849 29.16728 -0.696 0.487912
## cloacal_temp_C:week3 1.82229 0.51531 3.536 0.000591 ***
## cloacal_temp_C:week12 0.71647 0.98516 0.727 0.468586
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.613 on 112 degrees of freedom
## (1 observation deleted due to missingness)
## Multiple R-squared: 0.1888, Adjusted R-squared: 0.1526
## F-statistic: 5.213 on 5 and 112 DF, p-value: 0.0002434CEWL_cov_check_msmttemp <- lm(data = dat,
CEWL_g_m2h ~ msmt_temp_C*week)
summary(CEWL_cov_check_msmttemp)##
## Call:
## lm(formula = CEWL_g_m2h ~ msmt_temp_C * week, data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.2517 -1.5583 -0.3139 2.0494 10.7491
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 10.25244 3.27104 3.134 0.0022 **
## msmt_temp_C 0.02113 0.11155 0.189 0.8501
## week3 -60.68682 14.33616 -4.233 4.74e-05 ***
## week12 -3.76401 34.56181 -0.109 0.9135
## msmt_temp_C:week3 1.90326 0.46806 4.066 8.91e-05 ***
## msmt_temp_C:week12 0.16788 1.23504 0.136 0.8921
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.554 on 112 degrees of freedom
## (1 observation deleted due to missingness)
## Multiple R-squared: 0.2149, Adjusted R-squared: 0.1799
## F-statistic: 6.133 on 5 and 112 DF, p-value: 4.648e-05##
## Call:
## lm(formula = CEWL_g_m2h ~ msmt_VPD_kPa * week, data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.1892 -1.7946 -0.3724 2.0215 10.8386
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 11.13212 1.80429 6.170 1.12e-08 ***
## msmt_VPD_kPa -0.07911 0.52156 -0.152 0.880
## week3 -25.66036 5.31654 -4.827 4.43e-06 ***
## week12 0.33315 13.99194 0.024 0.981
## msmt_VPD_kPa:week3 6.37242 1.44871 4.399 2.49e-05 ***
## msmt_VPD_kPa:week12 0.20939 5.87420 0.036 0.972
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.504 on 112 degrees of freedom
## (1 observation deleted due to missingness)
## Multiple R-squared: 0.2368, Adjusted R-squared: 0.2027
## F-statistic: 6.95 on 5 and 112 DF, p-value: 1.093e-05## Anova Table (Type II tests)
##
## Response: CEWL_g_m2h
## Sum Sq Df F value Pr(>F)
## cloacal_temp_C 12.30 1 0.9422 0.333793
## week 169.15 2 6.4797 0.002173 **
## cloacal_temp_C:week 167.73 2 6.4253 0.002282 **
## Residuals 1461.83 112
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## Anova Table (Type II tests)
##
## Response: CEWL_g_m2h
## Sum Sq Df F value Pr(>F)
## msmt_temp_C 18.24 1 1.4441 0.2320199
## week 178.22 2 7.0546 0.0013014 **
## msmt_temp_C:week 208.88 2 8.2684 0.0004474 ***
## Residuals 1414.73 112
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## Anova Table (Type II tests)
##
## Response: CEWL_g_m2h
## Sum Sq Df F value Pr(>F)
## msmt_VPD_kPa 28.80 1 2.3452 0.1284889
## week 188.86 2 7.6898 0.0007423 ***
## msmt_VPD_kPa:week 237.73 2 9.6798 0.0001325 ***
## Residuals 1375.33 112
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1CEWL had a positive relationship with each of body temp, ambient temp, and ambient VPD, but only for the measurements taken in May (week 3).
Variance
Let’s also check the distribution/range of these covariates, because based on the figures, week 3 had a pretty small range of temperatures and VPDs.
dat %>%
group_by(week) %>%
summarise(ctemp_range = max(cloacal_temp_C, na.rm = T) -
min(cloacal_temp_C, na.rm = T),
mtemp_range = max(msmt_temp_C, na.rm = T) -
min(msmt_temp_C, na.rm = T),
mVPD_range = max(msmt_VPD_kPa, na.rm = T) -
min(msmt_VPD_kPa, na.rm = T))## # A tibble: 3 × 4
## week ctemp_range mtemp_range mVPD_range
## <fct> <dbl> <dbl> <dbl>
## 1 1 14 14.5 2.95
## 2 3 5 4.55 1.28
## 3 12 3 2.22 0.495Well… the ranges for week 3 are not small enough to be meaningless. So, we will have to keep in mind that temperature and VPD should be retained in CEWL models as a covariate.
Water Balance Metrics
What is the interrelation among water balance metrics CEWL, osmolality, hematocrit, and mass/SMI. I’ll use all measurements for this.
Models
Check a simple linear model between each set of paired measurements. I started with full models that included interactions with month and sex for each relationship, then dropped if the interaction was not significant. So, the models that still have interaction effects have them because they are important, and the models that do not have interaction effects had those effects dropped because they were not important.
## Anova Table (Type II tests)
##
## Response: CEWL_g_m2h
## Sum Sq Df F value Pr(>F)
## osmolality_mmol_kg_mean 17.59 1 1.3402 0.24950
## month 125.33 2 4.7747 0.01027 *
## osmolality_mmol_kg_mean:month 72.58 2 2.7652 0.06734 .
## Residuals 1443.69 110
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1##
## Call:
## lm(formula = CEWL_g_m2h ~ osmolality_mmol_kg_mean * month, data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.3779 -2.0003 0.0237 1.9817 10.8246
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.266624 6.658883 1.091 0.2775
## osmolality_mmol_kg_mean 0.010368 0.018227 0.569 0.5706
## monthMay 21.295848 10.315446 2.064 0.0413 *
## monthJuly -2.159012 23.011966 -0.094 0.9254
## osmolality_mmol_kg_mean:monthMay -0.062972 0.027584 -2.283 0.0244 *
## osmolality_mmol_kg_mean:monthJuly 0.008481 0.064845 0.131 0.8962
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.623 on 110 degrees of freedom
## (3 observations deleted due to missingness)
## Multiple R-squared: 0.1548, Adjusted R-squared: 0.1164
## F-statistic: 4.03 on 5 and 110 DF, p-value: 0.00214## Anova Table (Type II tests)
##
## Response: CEWL_g_m2h
## Sum Sq Df F value Pr(>F)
## hematocrit_percent 136.0 1 9.8621 0.002149 **
## Residuals 1572.1 114
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## Anova Table (Type II tests)
##
## Response: CEWL_g_m2h
## Sum Sq Df F value Pr(>F)
## mass_g 22.14 1 1.4432 0.2321
## Residuals 1779.91 116## Anova Table (Type II tests)
##
## Response: CEWL_g_m2h
## Sum Sq Df F value Pr(>F)
## SMI 58.54 1 3.8949 0.05081 .
## Residuals 1743.51 116
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# osmolality ~
car::Anova(lm(data = dat,
osmolality_mmol_kg_mean ~ hematocrit_percent), type = "II")## Anova Table (Type II tests)
##
## Response: osmolality_mmol_kg_mean
## Sum Sq Df F value Pr(>F)
## hematocrit_percent 9602 1 15.072 0.0001738 ***
## Residuals 72622 114
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## Anova Table (Type II tests)
##
## Response: osmolality_mmol_kg_mean
## Sum Sq Df F value Pr(>F)
## mass_g 0 1 3e-04 0.9869
## Residuals 82223 114## Anova Table (Type II tests)
##
## Response: osmolality_mmol_kg_mean
## Sum Sq Df F value Pr(>F)
## SMI 3656 1 5.3045 0.02308 *
## Residuals 78568 114
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## Anova Table (Type II tests)
##
## Response: hematocrit_percent
## Sum Sq Df F value Pr(>F)
## SMI 397.6 1 8.8942 0.003511 **
## sex_M_F 468.8 1 10.4855 0.001583 **
## SMI:sex_M_F 356.3 1 7.9694 0.005632 **
## Residuals 5007.3 112
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Notes on these relationships:
- CEWL ~ osml had significant interaction of osml* month (osml* mo* sex was not sig)
- CEWL ~ hct was sig but had no extra interactions (plot direct relationship only)
- CEWL ~ SMI was sig (mass was not; & SMI did not interact with mo or sex)
- osml & hct corr, but with no interactions of mo or sex
- osml also corr with SMI but not mass (& no interactions with mo or sex)
- ^all of the above based on anova SS, so assessing whether these things explain variation in response (CEWL/osml)
For MS:
Plasma osmolality explained a significant amount of the variation in CEWL (SLR: F(1,110) = 5.07, p = 0.03), with a marginally nonsignificant interaction with measurement time period (SLR: F(2,110) = 2.77, p = 0.07).
CEWL could also be explained by hematocrit (SLR: F(1,114) = 9.86, p = 0.002), or body condition (SLR: F(1,116) = 3.89, p = 0.05), but not lizard mass (SLR: F(1,116) = 1.44, p = 0.2).
Variation in plasma osmolality could be explained by hematocrit (SLR: F(1,114) = 15.07, p < 0.001) or body condition (SLR: F(1,114) = 5.30, p = 0.02), but not lizard mass (SLR: F(1,114) < 0.001, p = 1).
Hematocrit and body condition were also significantly related (SLR: F(1,112) = 9.03, p = 0.003), with an interaction of lizard sex (SLR: F(1,112) = 7.97, p = 0.005).
Plots
CEWL ~ Osmolality
ggplot(dat) +
aes(x = osmolality_mmol_kg_mean,
y = CEWL_g_m2h,
color = month,
shape = month
) +
geom_point(alpha = 0.5) +
stat_smooth(formula = y ~ x,
method = "lm",
se = F,
alpha = 0.8) +
scale_color_manual(values = n3, name = NULL
#name = "Measurement\nPeriod"
) +
scale_shape_manual(values = c(15, 16, 17), name = NULL) +
labs(x = bquote('Plasma Osmolality (mmol '*kg^-1*')'),
y = bquote('CEWL (g '*m^-2*' '*h^-1*')')) +
scale_x_continuous(limits = c(300, 450),
breaks = c(300, 350, 400, 450),
labels = c(300, 350, 400, 450)) +
savs_ggplot_theme -> CEWL_osml_fig
CEWL_osml_fig## Warning: Removed 3 rows containing non-finite values (`stat_smooth()`).## Warning: Removed 3 rows containing missing values (`geom_point()`).
# arrange to get legend centered
ggarrange(CEWL_osml_fig,
ncol = 1, nrow = 1,
common.legend = TRUE,
legend = "bottom"
) -> CEWL_osml_fig_arrange## Warning: Removed 3 rows containing non-finite values (`stat_smooth()`).
## Removed 3 rows containing missing values (`geom_point()`).## Warning: Removed 3 rows containing non-finite values (`stat_smooth()`).## Warning: Removed 3 rows containing missing values (`geom_point()`).CEWL ~ SMI
ggplot(dat) +
aes(x = SMI,
y = CEWL_g_m2h,
#color = tmt
#color = month_sex
) +
geom_point(size = 1,
alpha = 0.4) +
stat_smooth(formula = y ~ x,
method = "lm",
se = F,
size = 1.6,
alpha = 1,
color = n1) +
xlab("Body Condition (g')") +
ylab(bquote('CEWL (g/'*m^2*'h)')) +
savs_ggplot_theme -> CEWL_SMI_fig## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.## Warning: Removed 1 rows containing non-finite values (`stat_smooth()`).## Warning: Removed 1 rows containing missing values (`geom_point()`).
CEWL ~ Hematocrit
ggplot(dat) +
aes(x = hematocrit_percent,
y = CEWL_g_m2h,
#color = tmt
#color = month_sex
) +
geom_point(size = 1,
alpha = 0.4) +
stat_smooth(formula = y ~ x,
method = "lm",
se = F,
size = 1.6,
alpha = 1,
color = n1) +
xlab("Hematocrit (%)") +
ylab(bquote('CEWL (g/'*m^2*'h)')) +
#ylab("") +
#xlim(300, 450) +
#annotate(geom = "text", x = 430, y = 24, label = bquote(''*R^2*'=0.033'),
# size = 7) +
#ylim(0, 25) +
savs_ggplot_theme -> CEWL_hct_fig
CEWL_hct_fig## Warning: Removed 3 rows containing non-finite values (`stat_smooth()`).## Warning: Removed 3 rows containing missing values (`geom_point()`).
Osmol ~ SMI
ggplot(dat) +
aes(x = SMI,
y = osmolality_mmol_kg_mean,
#color = tmt
#color = month_sex
) +
geom_point(size = 1,
alpha = 0.4) +
stat_smooth(formula = y ~ x,
method = "lm",
se = F,
size = 1.6,
alpha = 1,
color = n1) +
xlab("Body Condition (g')") +
ylab("Osmolality (mmol/kg)") +
#ylab("") +
#xlim(300, 450) +
#annotate(geom = "text", x = 430, y = 24, label = bquote(''*R^2*'=0.033'),
# size = 7) +
#ylim(0, 25) +
savs_ggplot_theme -> osml_SMI_fig
osml_SMI_fig## Warning: Removed 3 rows containing non-finite values (`stat_smooth()`).## Warning: Removed 3 rows containing missing values (`geom_point()`).
Osmol ~ Hematocrit
ggplot(dat) +
aes(x = hematocrit_percent,
y = osmolality_mmol_kg_mean,
#color = tmt
#color = month_sex
) +
geom_point(size = 1,
alpha = 0.4) +
stat_smooth(formula = y ~ x,
method = "lm",
se = F,
size = 1.6,
alpha = 1,
color = n1) +
xlab("Hematocrit (%)") +
ylab("Osmolality (mmol/kg)") +
#ylab("") +
#xlim(300, 450) +
#annotate(geom = "text", x = 430, y = 24, label = bquote(''*R^2*'=0.033'),
# size = 7) +
#ylim(0, 25) +
scale_color_manual(values = n1) +
savs_ggplot_theme -> osml_hct_fig
osml_hct_fig## Warning: Removed 3 rows containing non-finite values (`stat_smooth()`).## Warning: Removed 3 rows containing missing values (`geom_point()`).
Hct ~ SMI
ggplot(dat) +
aes(x = SMI,
y = hematocrit_percent,
color = sex_M_F,
#color = month_sex
) +
geom_point(size = 1,
alpha = 0.4) +
stat_smooth(formula = y ~ x,
method = "lm",
se = F,
size = 1.6,
alpha = 1) +
xlab("Body Condition (g')") +
ylab("Hematocrit (%)") +
#ylab("") +
#xlim(300, 450) +
#annotate(geom = "text", x = 430, y = 24, label = bquote(''*R^2*'=0.033'),
# size = 7) +
#ylim(0, 25) +
scale_color_manual(values = n2) +
savs_ggplot_theme -> hct_SMI_fig
hct_SMI_fig## Warning: Removed 3 rows containing non-finite values (`stat_smooth()`).## Warning: Removed 3 rows containing missing values (`geom_point()`).
Supplemental Hydration Treatment Effects
We tried to give the lizards water to drink. Did they actually drink any? And did it affect measurable change in their physiology?
Calculate Deltas
First, calculate the change in CEWL, osmolality, and mass for each treatment individual that was recaptured between April and May.
# get the intitial values from April
starting_values <- dat_repeats %>%
# look only at lizards in hydration or sham tmt
dplyr::filter(tmt != "No Tmt") %>%
# get the values for April measurements only
dplyr::filter(month == "April") %>%
dplyr::select(individual_ID,
CEWL_April = CEWL_g_m2h,
osml_April = osmolality_mmol_kg_mean,
mass_April = mass_g)
# calculate deltas
dat_tmt_change <- dat_repeats %>%
# look only at lizards in hydration or sham tmt
dplyr::filter(tmt != "No Tmt") %>%
# add starting values
left_join(starting_values, by = 'individual_ID') %>%
# make sure we have complete data
dplyr::filter(complete.cases(mass_April, CEWL_April)) %>%
# calculate difference for each msmt versus in April
mutate(mass_change_from_April = mass_g - mass_April,
CEWL_change_from_April = CEWL_g_m2h - CEWL_April,
osml_change_from_April = osmolality_mmol_kg_mean - osml_April)
# subset of data to look at change April->May only
dat_tmt_change_May_only <- dat_tmt_change %>%
dplyr::filter(month == "May")Tmt Effect Check
What was mass change during treatment, when we either gave lizards water to drink or a sham treatment?
This quantifies the difference in mass gain/loss when we gave them water to drink versus sham treatment.
tmt_effect_mass <- dat %>%
dplyr::filter(tmt != "No Tmt") %>%
dplyr::filter(complete.cases(tmt_mass_change_g))
summary(tmt_effect_mass)## capture_date_time capture_date radio_collar_serial
## Min. :2021-04-23 03:39:00.00 Min. :2021-04-23 Length:61
## 1st Qu.:2021-04-23 07:54:00.00 1st Qu.:2021-04-23 Class :character
## Median :2021-04-24 03:49:00.00 Median :2021-04-24 Mode :character
## Mean :2021-04-28 02:55:24.21 Mean :2021-04-28
## 3rd Qu.:2021-05-07 04:31:00.00 3rd Qu.:2021-05-07
## Max. :2021-05-08 07:19:00.00 Max. :2021-05-08
## NA's :4
## PIT_tag_ID individual_ID sex_M_F mass_g SVL_mm
## Length:61 F-01 : 2 Female:29 Min. :26.00 Min. : 85
## Class :character F-05 : 2 Male :32 1st Qu.:33.00 1st Qu.: 95
## Mode :character F-06 : 2 Median :38.00 Median : 99
## F-11 : 2 Mean :39.38 Mean :101
## F-12 : 2 3rd Qu.:46.00 3rd Qu.:108
## F-13 : 2 Max. :56.00 Max. :122
## (Other):49
## tmt tmt_mass_change_g capture_to_msmt hematocrit_percent
## Water Tmt:29 Min. :-3.0000 Min. : 11.0 Min. :17.00
## Sham Tmt :32 1st Qu.:-1.0000 1st Qu.: 62.0 1st Qu.:29.75
## No Tmt : 0 Median : 0.0000 Median : 150.0 Median :36.00
## Mean :-0.2262 Mean : 403.5 Mean :34.30
## 3rd Qu.: 0.0000 3rd Qu.: 266.0 3rd Qu.:39.00
## Max. : 3.0000 Max. :1665.0 Max. :58.00
## NA's :4 NA's :1
## cloacal_temp_C CEWL_g_m2h msmt_temp_C msmt_RH_percent
## Min. :24.00 Min. : 0.650 Min. :19.07 Min. :12.20
## 1st Qu.:32.00 1st Qu.: 6.490 1st Qu.:28.57 1st Qu.:14.53
## Median :33.00 Median : 9.318 Median :30.50 Median :18.73
## Mean :32.44 Mean : 9.410 Mean :29.43 Mean :19.45
## 3rd Qu.:34.00 3rd Qu.:10.827 3rd Qu.:31.86 3rd Qu.:22.46
## Max. :36.00 Max. :21.673 Max. :33.55 Max. :35.83
##
## e_s_kPa msmt_VPD_kPa osmolality_mmol_kg_mean month week
## Min. :2.205 Min. :1.415 Min. :309.3 April:39 1 :39
## 1st Qu.:3.907 1st Qu.:3.014 1st Qu.:346.1 May :22 3 :22
## Median :4.365 Median :3.685 Median :364.2 July : 0 12: 0
## Mean :4.168 Mean :3.396 Mean :366.5
## 3rd Qu.:4.717 3rd Qu.:4.028 3rd Qu.:381.1
## Max. :5.188 Max. :4.319 Max. :436.5
## NA's :1
## week_num ultrasound_date gravid_Y_N number_eggs
## Min. :1.000 Min. :2021-04-24 Not Gravid:12 Min. :3.000
## 1st Qu.:1.000 1st Qu.:2021-04-24 Gravid :15 1st Qu.:3.500
## Median :1.000 Median :2021-04-24 NA's :34 Median :4.000
## Mean :1.721 Mean :2021-04-29 Mean :3.733
## 3rd Qu.:3.000 3rd Qu.:2021-05-08 3rd Qu.:4.000
## Max. :3.000 Max. :2021-05-08 Max. :4.000
## NA's :34 NA's :46
## largest_egg_diameter_cm largest_egg_circumference_cm stage
## Min. :0.570 Min. :1.710 small round: 1
## 1st Qu.:0.815 1st Qu.:2.320 large round: 9
## Median :0.990 Median :2.880 soft : 2
## Mean :1.073 Mean :3.045 soft oblong: 0
## 3rd Qu.:1.350 3rd Qu.:3.970 firm oblong: 2
## Max. :1.800 Max. :4.740 hard oblong: 1
## NA's :46 NA's :46 NA's :46
## dev_point SMI month_sex date_chunk_pretty
## Min. :1.000 Min. :27.04 April Female:19 Apr 23-25:39
## 1st Qu.:2.000 1st Qu.:32.99 April Male :20 May 7-9 :22
## Median :2.000 Median :34.76 May Female :10 Jul 14 : 0
## Mean :2.667 Mean :35.40 May Male :12
## 3rd Qu.:3.000 3rd Qu.:38.39 July Female : 0
## Max. :5.000 Max. :45.63 July Male : 0
## NA's :46
## percent_water_mass_change
## Min. :-6.8337
## 1st Qu.:-2.7027
## Median : 0.0000
## Mean :-0.5594
## 3rd Qu.: 0.0000
## Max. :10.0000
## # change in grams (actual amount of water "drank")
tmt_change <- lm(data = tmt_effect_mass, tmt_mass_change_g ~ tmt)
summary(tmt_change)##
## Call:
## lm(formula = tmt_mass_change_g ~ tmt, data = tmt_effect_mass)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.2344 -0.3690 -0.2344 0.6656 3.7656
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.3690 0.2090 1.765 0.082666 .
## tmtSham Tmt -1.1346 0.2886 -3.932 0.000224 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.125 on 59 degrees of freedom
## Multiple R-squared: 0.2076, Adjusted R-squared: 0.1942
## F-statistic: 15.46 on 1 and 59 DF, p-value: 0.0002241## tmt emmean SE df lower.CL upper.CL
## Water Tmt 0.369 0.209 59 -0.0492 0.787
## Sham Tmt -0.766 0.199 59 -1.1637 -0.368
##
## Confidence level used: 0.95# change in percent body mass (relative amount)
percent_change <- lm(data = tmt_effect_mass, percent_water_mass_change ~ tmt)
summary(percent_change)##
## Call:
## lm(formula = percent_water_mass_change ~ tmt, data = tmt_effect_mass)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.4465 -1.4239 -0.8996 1.8787 11.9094
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.9303 0.5617 1.656 0.102997
## tmtSham Tmt -2.8397 0.7756 -3.661 0.000538 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.025 on 59 degrees of freedom
## Multiple R-squared: 0.1852, Adjusted R-squared: 0.1713
## F-statistic: 13.41 on 1 and 59 DF, p-value: 0.0005376## tmt emmean SE df lower.CL upper.CL
## Water Tmt 0.93 0.562 59 -0.194 2.054
## Sham Tmt -1.91 0.535 59 -2.979 -0.839
##
## Confidence level used: 0.95So, it did result in a significant change in mass within the same-day water supplementation tmt (vs sham), but only ~1.4g difference, on average.
Change Over Time
Did hydric physiology change throughout these lizards’ active season? Based on hydration treatment? Sex?
Plot Change ~ Tmt
Mass
ggplot(dat_tmt_change) +
aes(x = month,
y = mass_change_from_April,
group = individual_ID,
color = tmt
) +
xlab("Measurement Month") +
ylab(expression(Delta ~ 'Lizard Mass (g)')) +
geom_line(alpha = 0.5) +
geom_point(size = 1.5) +
scale_color_manual(values = my_colors, name = "Sex") +
savs_ggplot_theme -> tmt_delta_mass
tmt_delta_mass
Osmolality
ggplot(dat_tmt_change) +
aes(x = month,
y = osml_change_from_April,
group = individual_ID,
color = tmt
) +
xlab("Measurement Month") +
ylab(expression(Delta ~ 'Osmolality (mmol '*kg^-1*')')) +
geom_line(alpha = 0.5) +
geom_point(size = 1.5) +
scale_color_manual(values = my_colors, name = "Sex") +
savs_ggplot_theme -> tmt_delta_osml
tmt_delta_osml## Warning: Removed 2 rows containing missing values (`geom_line()`).## Warning: Removed 2 rows containing missing values (`geom_point()`).
CEWL
ggplot(dat_tmt_change) +
aes(x = month,
y = CEWL_change_from_April,
group = individual_ID,
color = tmt
) +
xlab("Measurement Month") +
ylab(expression(Delta ~ 'CEWL (g '*m^-2*' '*h^-1*')')) +
geom_line(alpha = 0.5) +
geom_point(size = 1.5) +
scale_color_manual(values = my_colors, name = "Sex") +
savs_ggplot_theme -> tmt_delta_CEWL
tmt_delta_CEWL
It does not look like treatment led to any differences in the amount or direction of change in lizard mass, CEWL, or osmolality. Changes were also not consistent by sex.
Export Multi-Panel Fig
Save plots:
ggarrange(tmt_delta_CEWL,
tmt_delta_osml,
tmt_delta_mass,
ncol = 1, nrow = 3,
labels = c("A", "B", "C"),
common.legend = TRUE,
legend = "bottom"
) -> tmt_delta_multi_fig## Warning: Removed 2 rows containing missing values (`geom_line()`).## Warning: Removed 2 rows containing missing values (`geom_point()`).
Tmt Marginal Means
I want to know whether the average values for the lizards we measured changes over time based on supplemental hydration treatment. To test this, I can run a linear model, then use emmeans to get the marginal means and confidence intervals for each time point.
Do this only for lizards with repeat measurements, because those are the ones that received tmt.
Mass
## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## Missing cells for: monthJuly:tmtNo Tmt.
## Interpret type III hypotheses with care.## Type III Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 274.074 137.037 2 37.802 7.4633 0.001854 **
## tmt 224.282 112.141 2 32.029 6.1078 0.005661 **
## month:tmt 37.189 12.396 3 37.177 0.6751 0.572776
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## tmt = Water Tmt:
## month emmean SE df lower.CL upper.CL
## April 41.0 2.03 45.7 37.0 45.1
## May 46.1 2.26 56.7 41.6 50.6
## July 37.6 2.34 59.7 32.9 42.3
##
## tmt = Sham Tmt:
## month emmean SE df lower.CL upper.CL
## April 37.9 1.92 43.5 34.1 41.8
## May 39.7 1.99 47.8 35.7 43.7
## July 34.9 2.53 65.4 29.9 40.0
##
## tmt = No Tmt:
## month emmean SE df lower.CL upper.CL
## April 29.8 3.33 43.5 23.1 36.5
## May 31.8 3.33 43.5 25.1 38.5
## July nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95## tmt = Water Tmt:
## contrast estimate SE df t.ratio p.value
## April - May -5.05 1.99 37.6 -2.542 0.0395
## April - July 3.45 2.09 38.3 1.652 0.2369
## May - July 8.50 2.44 41.5 3.480 0.0033
##
## tmt = Sham Tmt:
## contrast estimate SE df t.ratio p.value
## April - May -1.76 1.65 35.7 -1.066 0.5412
## April - July 2.99 2.27 38.4 1.314 0.3964
## May - July 4.75 2.41 39.9 1.970 0.1328
##
## tmt = No Tmt:
## contrast estimate SE df t.ratio p.value
## April - May -2.00 2.71 35.0 -0.738 0.7428
## April - July nonEst NA NA NA NA
## May - July nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates## month = April:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt 3.12 2.79 44.7 1.116 0.5096
## Water Tmt - No Tmt 11.25 3.90 44.1 2.885 0.0163
## Sham Tmt - No Tmt 8.13 3.85 43.5 2.115 0.0984
##
## month = May:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt 6.40 3.01 52.9 2.127 0.0940
## Water Tmt - No Tmt 14.30 4.02 47.7 3.554 0.0024
## Sham Tmt - No Tmt 7.89 3.88 44.7 2.033 0.1159
##
## month = July:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt 2.65 3.45 63.4 0.769 0.7234
## Water Tmt - No Tmt nonEst NA NA NA NA
## Sham Tmt - No Tmt nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimatesBody Condition
## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## Missing cells for: monthJuly:tmtNo Tmt.
## Interpret type III hypotheses with care.## Type III Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 282.941 141.470 2 43.652 10.7878 0.0001559 ***
## tmt 42.122 21.061 2 33.543 1.6068 0.2155784
## month:tmt 45.205 15.068 3 41.791 1.1490 0.3405994
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## tmt = Water Tmt:
## month emmean SE df lower.CL upper.CL
## April 34.7 1.11 64.2 32.4 36.9
## May 39.6 1.34 65.9 37.0 42.3
## July 31.2 1.42 66.0 28.4 34.1
##
## tmt = Sham Tmt:
## month emmean SE df lower.CL upper.CL
## April 34.2 1.04 63.4 32.2 36.3
## May 36.2 1.11 64.8 33.9 38.4
## July 31.5 1.64 65.3 28.2 34.8
##
## tmt = No Tmt:
## month emmean SE df lower.CL upper.CL
## April 33.4 1.79 63.4 29.8 37.0
## May 34.4 1.79 63.4 30.8 38.0
## July nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95## tmt = Water Tmt:
## contrast estimate SE df t.ratio p.value
## April - May -4.99 1.62 43.0 -3.077 0.0100
## April - July 3.42 1.69 44.8 2.022 0.1186
## May - July 8.41 1.89 54.6 4.439 0.0001
##
## tmt = Sham Tmt:
## contrast estimate SE df t.ratio p.value
## April - May -1.94 1.38 37.3 -1.400 0.3514
## April - July 2.70 1.84 46.7 1.471 0.3141
## May - July 4.64 1.91 51.4 2.430 0.0481
##
## tmt = No Tmt:
## contrast estimate SE df t.ratio p.value
## April - May -1.03 2.29 35.5 -0.449 0.8953
## April - July nonEst NA NA NA NA
## May - July nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates## month = April:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt 0.429 1.52 63.9 0.282 0.9572
## Water Tmt - No Tmt 1.265 2.11 63.6 0.599 0.8211
## Sham Tmt - No Tmt 0.837 2.07 63.4 0.404 0.9142
##
## month = May:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt 3.479 1.74 65.6 1.997 0.1210
## Water Tmt - No Tmt 5.224 2.24 64.6 2.332 0.0584
## Sham Tmt - No Tmt 1.745 2.11 63.8 0.826 0.6881
##
## month = July:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt -0.291 2.17 65.8 -0.134 0.9901
## Water Tmt - No Tmt nonEst NA NA NA NA
## Sham Tmt - No Tmt nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimatesOsmolality
tmt_diffs_osml <- lmerTest::lmer(data = dat_repeats,
osmolality_mmol_kg_mean ~ month*tmt +
(1|individual_ID))## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## Missing cells for: monthJuly:tmtNo Tmt.
## Interpret type III hypotheses with care.## Type III Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 4511.6 2255.81 2 43.087 4.6155 0.01527 *
## tmt 343.7 171.85 2 33.790 0.3518 0.70600
## month:tmt 823.2 274.41 3 41.048 0.5614 0.64352
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## tmt = Water Tmt:
## month emmean SE df lower.CL upper.CL
## April 356 7.18 62.4 342 371
## May 379 8.29 63.8 363 396
## July 345 8.80 64.0 327 362
##
## tmt = Sham Tmt:
## month emmean SE df lower.CL upper.CL
## April 364 6.42 60.8 351 377
## May 375 6.90 62.5 361 389
## July 359 10.14 63.2 339 379
##
## tmt = No Tmt:
## month emmean SE df lower.CL upper.CL
## April 368 12.45 62.6 343 392
## May 376 11.13 60.8 354 399
## July nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95## tmt = Water Tmt:
## contrast estimate SE df t.ratio p.value
## April - May -23.01 10.04 39.7 -2.291 0.0687
## April - July 11.58 10.62 45.6 1.090 0.5249
## May - July 34.59 11.65 52.8 2.969 0.0122
##
## tmt = Sham Tmt:
## contrast estimate SE df t.ratio p.value
## April - May -10.87 8.45 35.5 -1.286 0.4122
## April - July 4.80 11.25 44.2 0.427 0.9048
## May - July 15.67 11.71 48.6 1.338 0.3812
##
## tmt = No Tmt:
## contrast estimate SE df t.ratio p.value
## April - May -8.74 15.05 37.3 -0.580 0.8314
## April - July nonEst NA NA NA NA
## May - July nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates## month = April:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt -7.53 9.64 61.8 -0.781 0.7159
## Water Tmt - No Tmt -11.11 14.38 62.5 -0.773 0.7208
## Sham Tmt - No Tmt -3.58 14.01 62.3 -0.256 0.9646
##
## month = May:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt 4.61 10.79 63.4 0.427 0.9045
## Water Tmt - No Tmt 3.16 13.88 62.2 0.227 0.9719
## Sham Tmt - No Tmt -1.45 13.09 61.4 -0.111 0.9933
##
## month = July:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt -14.31 13.42 63.8 -1.066 0.5384
## Water Tmt - No Tmt nonEst NA NA NA NA
## Sham Tmt - No Tmt nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimatesCEWL
## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## Missing cells for: monthJuly:tmtNo Tmt.
## Interpret type III hypotheses with care.## Type III Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 211.726 105.863 2 42.321 9.5390 0.0003805 ***
## tmt 39.478 19.739 2 34.064 1.7790 0.1841270
## month:tmt 43.448 14.483 3 40.735 1.3049 0.2858705
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# check covariates
tmt_diffs_CEWL_check_covaries <- lmerTest::lmer(data = dat_repeats,
CEWL_g_m2h ~ month*tmt
+ msmt_temp_C + msmt_VPD_kPa +
(1|individual_ID))## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient## Missing cells for: monthJuly:tmtNo Tmt.
## Interpret type III hypotheses with care.## Type III Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 158.396 79.198 2 46.720 6.8694 0.002425 **
## tmt 35.565 17.783 2 33.062 1.5441 0.228516
## msmt_temp_C 5.959 5.959 1 61.410 0.5174 0.474666
## msmt_VPD_kPa 7.641 7.641 1 61.763 0.6635 0.418456
## month:tmt 42.580 14.193 3 39.151 1.2321 0.310989
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## tmt = Water Tmt:
## month emmean SE df lower.CL upper.CL
## April 9.32 1.09 60.8 7.13 11.50
## May 6.11 1.30 64.6 3.51 8.71
## July 11.60 1.38 64.9 8.84 14.35
##
## tmt = Sham Tmt:
## month emmean SE df lower.CL upper.CL
## April 10.04 1.02 59.2 8.00 12.08
## May 9.30 1.09 61.9 7.12 11.48
## July 12.98 1.58 64.1 9.83 16.13
##
## tmt = No Tmt:
## month emmean SE df lower.CL upper.CL
## April 13.43 1.97 62.3 9.50 17.37
## May 7.87 1.77 59.2 4.34 11.41
## July nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95## tmt = Water Tmt:
## contrast estimate SE df t.ratio p.value
## April - May 3.208 1.51 40.6 2.131 0.0960
## April - July -2.278 1.57 42.2 -1.447 0.3265
## May - July -5.487 1.78 50.5 -3.075 0.0094
##
## tmt = Sham Tmt:
## contrast estimate SE df t.ratio p.value
## April - May 0.741 1.28 35.9 0.581 0.8312
## April - July -2.942 1.71 43.4 -1.719 0.2096
## May - July -3.683 1.79 47.3 -2.060 0.1093
##
## tmt = No Tmt:
## contrast estimate SE df t.ratio p.value
## April - May 5.558 2.28 37.9 2.439 0.0500
## April - July nonEst NA NA NA NA
## May - July nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates## month = April:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt -0.721 1.50 60.1 -0.482 0.8800
## Water Tmt - No Tmt -4.112 2.25 62.0 -1.826 0.1696
## Sham Tmt - No Tmt -3.391 2.22 61.7 -1.529 0.2845
##
## month = May:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt -3.188 1.70 63.8 -1.876 0.1538
## Water Tmt - No Tmt -1.763 2.20 61.7 -0.803 0.7027
## Sham Tmt - No Tmt 1.426 2.08 60.0 0.686 0.7723
##
## month = July:
## contrast estimate SE df t.ratio p.value
## Water Tmt - Sham Tmt -1.385 2.10 64.8 -0.661 0.7869
## Water Tmt - No Tmt nonEst NA NA NA NA
## Sham Tmt - No Tmt nonEst NA NA NA NA
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimatesHydric physiology was not significantly different between treatments at any time.
Was change overall (across tmts) different from zero?
##
## One Sample t-test
##
## data: dat_tmt_change_May_only$mass_change_from_April
## t = 3.828, df = 20, p-value = 0.001051
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## 1.259063 4.274270
## sample estimates:
## mean of x
## 2.766667##
## One Sample t-test
##
## data: dat_tmt_change_May_only$CEWL_change_from_April
## t = -1.8247, df = 20, p-value = 0.08303
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## -4.805881 0.321135
## sample estimates:
## mean of x
## -2.242373##
## One Sample t-test
##
## data: dat_tmt_change_May_only$osml_change_from_April
## t = 1.8095, df = 20, p-value = 0.08542
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## -2.027158 28.566840
## sample estimates:
## mean of x
## 13.26984On average, lizards who were measured in both April and May increased mass, but did not experience a change in osmolality or CEWL.
Plot Delta ~ Time*Sex
Mass
ggplot(dat_tmt_change) +
aes(x = month,
y = mass_change_from_April,
group = individual_ID,
color = sex_M_F,
shape = sex_M_F
) +
xlab("Measurement Month") +
ylab(expression(Delta ~ 'Lizard Mass (g)')) +
geom_line(alpha = 0.5) +
geom_point(size = 1.5) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes, name = "Sex") +
savs_ggplot_theme -> delta_mass_sex
delta_mass_sex
Osmolality
ggplot(dat_tmt_change) +
aes(x = month,
y = osml_change_from_April,
group = individual_ID,
color = sex_M_F,
shape = sex_M_F
) +
xlab("Measurement Month") +
ylab(expression(Delta ~ 'Osmolality (mmol '*kg^-1*')')) +
geom_line(alpha = 0.5) +
geom_point(size = 1.5) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes, name = "Sex") +
savs_ggplot_theme -> delta_osml_sex
delta_osml_sex## Warning: Removed 2 rows containing missing values (`geom_line()`).## Warning: Removed 2 rows containing missing values (`geom_point()`).
CEWL
ggplot(dat_tmt_change) +
aes(x = month,
y = CEWL_change_from_April,
group = individual_ID,
color = sex_M_F,
shape = sex_M_F
) +
xlab("Measurement Month") +
ylab(expression(Delta ~ 'CEWL (g '*m^-2*' '*h^-1*')')) +
geom_line(alpha = 0.5) +
geom_point(size = 1.5) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes, name = "Sex") +
savs_ggplot_theme -> delta_CEWL_sex
delta_CEWL_sex
Export Multi-Panel Fig
Save plots:
ggarrange(delta_CEWL_sex,
delta_osml_sex,
delta_mass_sex,
ncol = 3, nrow = 1,
#ncol = 1, nrow = 3,
labels = c("A", "B", "C"),
common.legend = TRUE,
legend = "bottom"
) -> delta_multi_fig_sex## Warning: Removed 2 rows containing missing values (`geom_line()`).## Warning: Removed 2 rows containing missing values (`geom_point()`).Time*Sex Marginal Means
I want to know whether the average values for the lizards we measured changes over time. To test this, I can run a linear model, then use emmeans to get the marginal means and confidence intervals for each time point.
Do this only for lizards with repeat measurements.
info for using emmeans with interactions: https://cran.r-project.org/web/packages/emmeans/vignettes/interactions.html
Check n’s:
## `summarise()` has grouped output by 'month'. You can override using the
## `.groups` argument.## # A tibble: 6 × 3
## # Groups: month [3]
## month sex_M_F n
## <fct> <fct> <int>
## 1 April Female 11
## 2 April Male 22
## 3 May Female 10
## 4 May Male 17
## 5 July Female 2
## 6 July Male 12Mass
time_sex_diffs_mass <- lmerTest::lmer(data = dat_repeats, mass_g ~ month*sex_M_F +
(1|individual_ID))
anova(time_sex_diffs_mass, type = "3", ddf = "Kenward-Roger")## Type III Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 314.003 157.002 2 38.639 9.6139 0.0004086 ***
## sex_M_F 0.032 0.032 1 38.343 0.0020 0.9647658
## month:sex_M_F 92.617 46.309 2 38.639 2.8357 0.0709404 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# mean values
emmeans_mass <- data.frame(emmeans(time_sex_diffs_mass,
pairwise ~ month | sex_M_F)$emmean) %>%
mutate(measurement = "Mass")
# stat diffs by month
time_emmeans_mass <- data.frame(emmeans(time_sex_diffs_mass,
pairwise ~ month | sex_M_F)$contrasts) %>%
mutate(measurement = "Mass")
# stat diffs by sex
sex_emmeans_mass <- data.frame(emmeans(time_sex_diffs_mass,
pairwise ~ sex_M_F | month)$contrasts) %>%
mutate(measurement = "Mass")Body Condition
time_sex_diffs_SMI <- lmerTest::lmer(data = dat_repeats, SMI ~ month*sex_M_F +
(1|individual_ID))
anova(time_sex_diffs_SMI, type = "3", ddf = "Kenward-Roger")## Type III Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 327.19 163.597 2 45.863 15.2615 8.308e-06 ***
## sex_M_F 0.63 0.628 1 46.566 0.0587 0.809671
## month:sex_M_F 134.56 67.280 2 45.863 6.2764 0.003896 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# mean values
emmeans_SMI <- data.frame(emmeans(time_sex_diffs_SMI,
pairwise ~ month | sex_M_F)$emmean) %>%
mutate(measurement = "SMI")
# stat diffs by month
time_emmeans_SMI <- data.frame(emmeans(time_sex_diffs_SMI,
pairwise ~ month | sex_M_F)$contrasts) %>%
mutate(measurement = "SMI")
# stat diffs by sex
sex_emmeans_SMI <- data.frame(emmeans(time_sex_diffs_SMI,
pairwise ~ sex_M_F | month)$contrasts) %>%
mutate(measurement = "SMI")Hematocrit
time_sex_diffs_hct <- lmerTest::lmer(data = dat_repeats,
hematocrit_percent ~ month*sex_M_F +
(1|individual_ID))
anova(time_sex_diffs_hct, type = "2", ddf = "Kenward-Roger")## Type II Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 817.60 408.80 2 44.421 15.8405 6.412e-06 ***
## sex_M_F 169.86 169.86 1 31.941 6.5855 0.01518 *
## month:sex_M_F 21.74 10.87 2 45.442 0.4211 0.65885
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# mean values
emmeans_hct <- data.frame(emmeans(time_sex_diffs_hct,
pairwise ~ month | sex_M_F)$emmean) %>%
mutate(measurement = "Hematocrit")
# stat diffs by month
time_emmeans_hct <- data.frame(emmeans(time_sex_diffs_hct,
pairwise ~ month | sex_M_F)$contrasts) %>%
mutate(measurement = "Hematocrit")
# stat diffs by sex
sex_emmeans_hct <- data.frame(emmeans(time_sex_diffs_hct,
pairwise ~ sex_M_F | month)$contrasts) %>%
mutate(measurement = "Hematocrit")Osmolality
time_sex_diffs_osml <- lmerTest::lmer(data = dat_repeats,
osmolality_mmol_kg_mean ~ month*sex_M_F +
(1|individual_ID))## boundary (singular) fit: see help('isSingular')## Type II Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 8073.9 4037.0 2 48.038 9.6102 0.0003086 ***
## sex_M_F 11642.2 11642.2 1 31.800 27.7401 9.302e-06 ***
## month:sex_M_F 1301.2 650.6 2 49.692 1.5476 0.2228392
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# mean values
emmeans_osml <- data.frame(emmeans(time_sex_diffs_osml,
pairwise ~ month | sex_M_F)$emmean) %>%
mutate(measurement = "Osmolality")
# stat diffs by month
time_emmeans_osml <- data.frame(emmeans(time_sex_diffs_osml,
pairwise ~ month | sex_M_F)$contrasts) %>%
mutate(measurement = "Osmolality")
# stat diffs by sex
sex_emmeans_osml <- data.frame(emmeans(time_sex_diffs_osml,
pairwise ~ sex_M_F | month)$contrasts) %>%
mutate(measurement = "Osmolality")CEWL
time_sex_diffs_CEWL <- lmerTest::lmer(data = dat_repeats,
CEWL_g_m2h ~ month*sex_M_F +
(1|individual_ID))
anova(time_sex_diffs_CEWL, type = "2", ddf = "Kenward-Roger")## Type II Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 180.393 90.197 2 44.819 8.0478 0.001032 **
## sex_M_F 103.055 103.055 1 32.419 9.1999 0.004737 **
## month:sex_M_F 18.921 9.461 2 46.365 0.8438 0.436564
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# check covariates
time_sex_diffs_CEWL_check_covaries <- lmerTest::lmer(data = dat_repeats,
CEWL_g_m2h ~ month*sex_M_F
+ msmt_temp_C + msmt_VPD_kPa +
(1|individual_ID))
anova(time_sex_diffs_CEWL_check_covaries, type = "2", ddf = "Kenward-Roger") # not important at all## Type II Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 101.714 50.857 2 50.265 4.3866 0.01754 *
## sex_M_F 89.529 89.529 1 34.970 7.7316 0.00868 **
## msmt_temp_C 0.133 0.133 1 64.720 0.0115 0.91496
## msmt_VPD_kPa 0.090 0.090 1 64.849 0.0078 0.93004
## month:sex_M_F 17.518 8.759 2 46.384 0.7555 0.47547
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## [1] 655.5556## [1] 385.7143# mean values
emmeans_CEWL <- data.frame(emmeans(time_sex_diffs_CEWL,
pairwise ~ month | sex_M_F)$emmean) %>%
mutate(measurement = "CEWL")
# stat diffs by month
time_emmeans_CEWL <- data.frame(emmeans(time_sex_diffs_CEWL,
pairwise ~ month | sex_M_F)$contrasts) %>%
mutate(measurement = "CEWL")
# stat diffs by sex
sex_emmeans_CEWL <- data.frame(emmeans(time_sex_diffs_CEWL,
pairwise ~ sex_M_F | month)$contrasts) %>%
mutate(measurement = "CEWL")Check Assumptions
Merge dfs
emmeans_all <- emmeans_CEWL %>%
rbind(emmeans_osml) %>%
rbind(emmeans_mass) %>%
rbind(emmeans_SMI) %>%
mutate(month_sex = factor(paste(month, sex_M_F),
levels = c("April Female", "April Male",
"May Female", "May Male",
"July Female", "July Male")))
# all pairwise stats together
emmean_time_stats <- time_emmeans_CEWL %>%
rbind(time_emmeans_osml) %>%
rbind(time_emmeans_mass) %>%
rbind(time_emmeans_SMI)
emmean_sex_stats <- sex_emmeans_CEWL %>%
rbind(sex_emmeans_osml) %>%
rbind(sex_emmeans_mass) %>%
rbind(sex_emmeans_SMI)Plot MMeans ~ Time*Sex
These are plots of the marginal means, only for the dataset of lizards with repeat measurements.
Mass
ggplot() +
geom_jitter(data = dat_repeats,
aes(x = month_sex,
y = mass_g,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = emmeans_all[emmeans_all$measurement == "Mass",],
aes(x = month_sex,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = emmeans_all[emmeans_all$measurement == "Mass",],
aes(x = month_sex,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab("Lizard Mass (g)") +
# sex differences = NONE
# time differences for MALES
annotate(geom = "text", x = 2, y = 61, label = "A", size = 3) +
annotate(geom = "text", x = 4, y = 61, label = "A", size = 3) +
annotate(geom = "text", x = 6, y = 61, label = "A", size = 3) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 61, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 61, label = "b", size = 3) +
annotate(geom = "text", x = 5, y = 61, label = "c", size = 3) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_fill_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes_open, name = "Sex") +
scale_x_discrete(labels = my_labels) +
ylim(17, 61) +
savs_ggplot_theme -> mass_emmean_plot
mass_emmean_plot
Body Condition
ggplot() +
geom_jitter(data = dat_repeats,
aes(x = month_sex,
y = SMI,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = emmeans_all[emmeans_all$measurement == "SMI",],
aes(x = month_sex,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = emmeans_all[emmeans_all$measurement == "SMI",],
aes(x = month_sex,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab("Body Condition (g')") +
# sex differences
annotate(geom = "text", x = 3.5, y = 45, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 2, y = 54, label = "AB", size = 3) +
annotate(geom = "text", x = 4, y = 54, label = "A", size = 3) +
annotate(geom = "text", x = 6, y = 54, label = "B", size = 3) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 54, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 54, label = "b", size = 3) +
annotate(geom = "text", x = 5, y = 54, label = "a", size = 3) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_fill_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes_open, name = "Sex") +
scale_x_discrete(labels = my_labels) +
ylim(20, 55) +
savs_ggplot_theme -> SMI_emmean_plot
SMI_emmean_plot
Osmolality
ggplot() +
geom_jitter(data = dat_repeats,
aes(x = month_sex,
y = osmolality_mmol_kg_mean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = emmeans_all[emmeans_all$measurement == "Osmolality",],
aes(x = month_sex,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = emmeans_all[emmeans_all$measurement == "Osmolality",],
aes(x = month_sex,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab(bquote('Osmolality (mmol '*kg^-1*')')) +
# sex differences
annotate(geom = "text", x = 1.5, y = 416, label = "*", size = 6) +
annotate(geom = "text", x = 3.5, y = 416, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 2, y = 430, label = "A", size = 3) +
annotate(geom = "text", x = 4, y = 430, label = "A", size = 3) +
annotate(geom = "text", x = 6, y = 430, label = "B", size = 3) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 430, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 430, label = "a", size = 3) +
annotate(geom = "text", x = 5, y = 430, label = "a", size = 3) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_fill_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes_open, name = "Sex") +
scale_x_discrete(labels = my_labels) +
scale_y_continuous(limits = c(320, 440), breaks = c(320, 360, 400, 440)) +
savs_ggplot_theme -> osml_emmean_plot
osml_emmean_plot## Warning: Removed 2 rows containing missing values (`geom_point()`).
CEWL
ggplot() +
geom_jitter(data = dat_repeats,
aes(x = month_sex,
y = CEWL_g_m2h,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = emmeans_all[emmeans_all$measurement == "CEWL",],
aes(x = month_sex,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = emmeans_all[emmeans_all$measurement == "CEWL",],
aes(x = month_sex,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
# sex differences
annotate(geom = "text", x = 3.5, y = 16, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 2, y = 19, label = "A", size = 3) +
annotate(geom = "text", x = 4, y = 19, label = "B", size = 3) +
annotate(geom = "text", x = 6, y = 19, label = "A", size = 3) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 19, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 19, label = "a", size = 3) +
annotate(geom = "text", x = 5, y = 19, label = "a", size = 3) +
xlab("") +
ylab(bquote('CEWL (g '*m^-2*' '*h^-1*')')) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_fill_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes_open, name = "Sex") +
scale_x_discrete(labels = my_labels) +
savs_ggplot_theme -> CEWL_emmean_plot
CEWL_emmean_plot## Warning: Removed 1 rows containing missing values (`geom_point()`).
This is only for lizards that were recaptured and measured more than once, so I think these changes may reflect true seasonality in their water balance.
Export Repeats-Only Multi-Panel Fig
Put in one plot together:
ggarrange(CEWL_emmean_plot,
osml_emmean_plot,
mass_emmean_plot,
SMI_emmean_plot,
ncol = 4, nrow = 1,
# ncol = 1, nrow = 3,
labels = c("A", "B", "C", "D"),
common.legend = TRUE,
legend = "bottom"
) -> szn_change_multi_fig## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Removed 1 rows containing missing values (`geom_point()`).## Warning: Removed 2 rows containing missing values (`geom_point()`).FULL Time*Sex Marginal Means
I’m going to do literally the same exact thing as the stats and figures in the section above, but now for all of my data, not just for the data on lizards with repeat measurements.
Check n’s:
## `summarise()` has grouped output by 'month'. You can override using the
## `.groups` argument.## # A tibble: 6 × 3
## # Groups: month [3]
## month sex_M_F n
## <fct> <fct> <int>
## 1 April Female 26
## 2 April Male 41
## 3 May Female 15
## 4 May Male 22
## 5 July Female 3
## 6 July Male 12Mass
full_time_sex_diffs_mass <- lmerTest::lmer(data = dat, mass_g ~ month*sex_M_F +
(1|individual_ID))
full_time_sex_diffs_mass_aov <- data.frame(anova(full_time_sex_diffs_mass,
type = "2", ddf = "Kenward-Roger")) %>%
mutate(explanatory = factor(row.names(.)))
# mean values
full_emmeans_mass <- data.frame(emmeans(full_time_sex_diffs_mass,
pairwise ~ month | sex_M_F)$emmean) %>%
mutate(measurement = "Mass")
# stat diffs by month
full_time_emmeans_mass <- data.frame(emmeans(full_time_sex_diffs_mass,
pairwise ~ month | sex_M_F)$contrasts) %>%
mutate(measurement = "Mass")
# stat diffs by sex
full_sex_emmeans_mass <- data.frame(emmeans(full_time_sex_diffs_mass,
pairwise ~ sex_M_F | month)$contrasts) %>%
mutate(measurement = "Mass")Body Condition
full_time_sex_diffs_SMI <- lmerTest::lmer(data = dat, SMI ~ month*sex_M_F +
(1|individual_ID))
full_time_sex_diffs_SMI_aov <- data.frame(anova(full_time_sex_diffs_SMI,
type = "3", ddf = "Kenward-Roger")) %>%
mutate(explanatory = factor(row.names(.)))
# double check that type 2 vs 3 anovas are essentially the same results
full_time_sex_diffs_SMI_aov2 <- data.frame(anova(full_time_sex_diffs_SMI,
type = "2", ddf = "Kenward-Roger")) %>%
mutate(explanatory = factor(row.names(.)))
full_time_sex_diffs_SMI_aov## Sum.Sq Mean.Sq NumDF DenDF F.value Pr..F.
## month 400.48274 200.24137 2 79.56813 18.583496 2.384786e-07
## sex_M_F 11.62646 11.62646 1 102.97466 1.079132 3.013248e-01
## month:sex_M_F 138.38877 69.19439 2 79.56813 6.421618 2.597665e-03
## explanatory
## month month
## sex_M_F sex_M_F
## month:sex_M_F month:sex_M_F## Sum.Sq Mean.Sq NumDF DenDF F.value Pr..F.
## month 348.170763 174.085381 2 71.15328 16.157819 1.639477e-06
## sex_M_F 8.594231 8.594231 1 74.72959 0.797690 3.746521e-01
## month:sex_M_F 138.388775 69.194387 2 79.56813 6.421618 2.597665e-03
## explanatory
## month month
## sex_M_F sex_M_F
## month:sex_M_F month:sex_M_F# mean values
full_emmeans_SMI <- data.frame(emmeans(full_time_sex_diffs_SMI,
pairwise ~ month | sex_M_F)$emmean) %>%
mutate(measurement = "SMI")
# stat diffs by month
full_time_emmeans_SMI <- data.frame(emmeans(full_time_sex_diffs_SMI,
pairwise ~ month | sex_M_F)$contrasts) %>%
mutate(measurement = "SMI")
# stat diffs by sex
full_sex_emmeans_SMI <- data.frame(emmeans(full_time_sex_diffs_SMI,
pairwise ~ sex_M_F | month)$contrasts) %>%
mutate(measurement = "SMI")Osmolality
full_time_sex_diffs_osml <- lmerTest::lmer(data = dat,
osmolality_mmol_kg_mean ~ month*sex_M_F +
(1|individual_ID))
full_time_sex_diffs_osml_aov <- data.frame(anova(full_time_sex_diffs_osml,
type = "3", ddf = "Kenward-Roger")) %>%
mutate(explanatory = factor(row.names(.)))
full_time_sex_diffs_osml_aov2 <- data.frame(anova(full_time_sex_diffs_osml,
type = "2", ddf = "Kenward-Roger")) %>%
mutate(explanatory = factor(row.names(.)))
full_time_sex_diffs_osml_aov## Sum.Sq Mean.Sq NumDF DenDF F.value Pr..F.
## month 6748.7057 3374.3528 2 91.86749 5.8321054 0.004125822
## sex_M_F 1441.1831 1441.1831 1 88.14015 2.4917205 0.118028949
## month:sex_M_F 770.2048 385.1024 2 91.86749 0.6655966 0.516427511
## explanatory
## month month
## sex_M_F sex_M_F
## month:sex_M_F month:sex_M_F## Sum.Sq Mean.Sq NumDF DenDF F.value Pr..F.
## month 9627.0309 4813.5154 2 87.00381 8.3213727 0.0004933249
## sex_M_F 5545.4928 5545.4928 1 69.65080 9.5878298 0.0028217467
## month:sex_M_F 770.2048 385.1024 2 91.86749 0.6655966 0.5164275115
## explanatory
## month month
## sex_M_F sex_M_F
## month:sex_M_F month:sex_M_F# mean values
full_emmeans_osml <- data.frame(emmeans(full_time_sex_diffs_osml,
pairwise ~ month | sex_M_F)$emmean) %>%
mutate(measurement = "Osmolality")
# stat diffs by month
full_time_emmeans_osml <- data.frame(emmeans(full_time_sex_diffs_osml,
pairwise ~ month | sex_M_F)$contrasts) %>%
mutate(measurement = "Osmolality")
# stat diffs by sex
full_sex_emmeans_osml <- data.frame(emmeans(full_time_sex_diffs_osml,
pairwise ~ sex_M_F | month)$contrasts) %>%
mutate(measurement = "Osmolality")CEWL
full_time_sex_diffs_CEWL <- lmerTest::lmer(data = dat,
CEWL_g_m2h ~ month*sex_M_F +
(1|individual_ID))
full_time_sex_diffs_CEWL_aov <- data.frame(anova(full_time_sex_diffs_CEWL,
type = "3", ddf = "Kenward-Roger")) %>%
mutate(explanatory = factor(row.names(.)))
full_time_sex_diffs_CEWL_aov2 <- data.frame(anova(full_time_sex_diffs_CEWL,
type = "2", ddf = "Kenward-Roger")) %>%
mutate(explanatory = factor(row.names(.)))
full_time_sex_diffs_CEWL_aov## Sum.Sq Mean.Sq NumDF DenDF F.value Pr..F.
## month 119.80851 59.90425 2 85.07588 5.625373 0.005075837
## sex_M_F 34.27453 34.27453 1 98.76587 3.219121 0.075842334
## month:sex_M_F 46.88928 23.44464 2 85.07588 2.201594 0.116890977
## explanatory
## month month
## sex_M_F sex_M_F
## month:sex_M_F month:sex_M_F## Sum.Sq Mean.Sq NumDF DenDF F.value Pr..F.
## month 181.78030 90.89015 2 77.41997 8.536304 0.0004465354
## sex_M_F 55.19526 55.19526 1 73.45087 5.184032 0.0257090750
## month:sex_M_F 46.88928 23.44464 2 85.07588 2.201594 0.1168909770
## explanatory
## month month
## sex_M_F sex_M_F
## month:sex_M_F month:sex_M_F# check covariates
full_time_sex_diffs_CEWL_check_covaries <- lmerTest::lmer(data = dat,
CEWL_g_m2h ~ month*sex_M_F
+ msmt_temp_C + msmt_VPD_kPa +
(1|individual_ID))
anova(full_time_sex_diffs_CEWL_check_covaries,
type = "2", ddf = "Kenward-Roger")## Type II Analysis of Variance Table with Kenward-Roger's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## month 185.755 92.878 2 89.095 8.5769 0.0003921 ***
## sex_M_F 46.450 46.450 1 74.944 4.2908 0.0417603 *
## msmt_temp_C 6.098 6.098 1 109.982 0.5633 0.4545481
## msmt_VPD_kPa 9.009 9.009 1 109.977 0.8322 0.3636216
## month:sex_M_F 32.730 16.365 2 86.160 1.5113 0.2264115
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# mean values
full_emmeans_CEWL <- data.frame(emmeans(full_time_sex_diffs_CEWL,
pairwise ~ month | sex_M_F)$emmean) %>%
mutate(measurement = "CEWL")
# stat diffs by month
full_time_emmeans_CEWL <- data.frame(emmeans(full_time_sex_diffs_CEWL,
pairwise ~ month | sex_M_F)$contrasts) %>%
mutate(measurement = "CEWL")
# stat diffs by sex
full_sex_emmeans_CEWL <- data.frame(emmeans(full_time_sex_diffs_CEWL,
pairwise ~ sex_M_F | month)$contrasts) %>%
mutate(measurement = "CEWL")Check Assumptions
Merge dfs
full_emmeans_all <- full_emmeans_CEWL %>%
rbind(full_emmeans_osml) %>%
rbind(full_emmeans_mass) %>%
rbind(full_emmeans_SMI) %>%
mutate(month_sex = factor(paste(month, sex_M_F),
levels = c("April Female", "April Male",
"May Female", "May Male",
"July Female", "July Male")))
full_emmeans_F <- full_emmeans_all %>%
dplyr::filter(sex_M_F == "Female")
full_emmeans_M <- full_emmeans_all %>%
dplyr::filter(sex_M_F == "Male")
# all pairwise stats together
full_emmean_time_stats <- full_time_emmeans_CEWL %>%
rbind(full_time_emmeans_osml) %>%
rbind(full_time_emmeans_mass) %>%
rbind(full_time_emmeans_SMI)
full_emmean_sex_stats <- full_sex_emmeans_CEWL %>%
rbind(full_sex_emmeans_osml) %>%
rbind(full_sex_emmeans_mass) %>%
rbind(full_sex_emmeans_SMI)
# anova stats - type 2 is better for looking at interactions
full_time_sex_diffs_all <- (full_time_sex_diffs_CEWL_aov2) %>%
rbind((full_time_sex_diffs_osml_aov2)) %>%
#rbind((full_time_sex_diffs_mass_aov2)) %>%
rbind((full_time_sex_diffs_SMI_aov2)) %>%
mutate(response = factor(c(rep("CEWL", 3), rep("Plasma Osmolality", 3),
#rep("Mass", 3),
rep("Body Condition", 3))))
summary(full_time_sex_diffs_all)## Sum.Sq Mean.Sq NumDF DenDF
## Min. : 8.594 Min. : 8.594 Min. :1.000 Min. :69.65
## 1st Qu.: 55.195 1st Qu.: 55.195 1st Qu.:1.000 1st Qu.:73.45
## Median : 181.780 Median : 90.890 Median :2.000 Median :77.42
## Mean :1857.972 Mean :1240.613 Mean :1.667 Mean :78.88
## 3rd Qu.: 770.205 3rd Qu.: 385.102 3rd Qu.:2.000 3rd Qu.:85.08
## Max. :9627.031 Max. :5545.493 Max. :2.000 Max. :91.87
## F.value Pr..F. explanatory response
## Min. : 0.6656 Min. :0.0000016 month :3 Body Condition :3
## 1st Qu.: 2.2016 1st Qu.:0.0004933 month:sex_M_F:3 CEWL :3
## Median : 6.4216 Median :0.0028217 sex_M_F :3 Plasma Osmolality:3
## Mean : 6.4304 Mean :0.1155601
## 3rd Qu.: 8.5363 3rd Qu.:0.1168910
## Max. :16.1578 Max. :0.5164275Export Model Results
full_time_sex_diffs_4exp <- full_time_sex_diffs_all %>%
mutate(partial_sum_of_squares = round(Sum.Sq, 0),
F_statistic = round(F.value, 2),
DenDF = round(DenDF, 0),
F_full = paste(F_statistic, " (", NumDF, ",", DenDF, ")", sep = ""),
p_value = round(Pr..F., 4)) %>%
dplyr::select(response, explanatory,
partial_sum_of_squares,
#NumDF, DenDF,
#F_statistic,
F_full,
p_value)
write.csv(full_time_sex_diffs_4exp,
"./results_statistics/hydric_Fstats_by_monthsex.csv")FULL Plot MMeans Time*Sex
Mass
ggplot() +
geom_jitter(data = dat,
aes(x = month_sex,
y = mass_g,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_all[full_emmeans_all$measurement == "Mass",],
aes(x = month_sex,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_all[full_emmeans_all$measurement == "Mass",],
aes(x = month_sex,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab("Lizard Mass (g)") +
# sex differences
annotate(geom = "text", x = 5.5, y = 54, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 2, y = 61, label = "AB", size = 3) +
annotate(geom = "text", x = 4, y = 61, label = "A", size = 3) +
annotate(geom = "text", x = 6, y = 61, label = "B", size = 3) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 61, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 61, label = "b", size = 3) +
annotate(geom = "text", x = 5, y = 61, label = "a", size = 3) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_fill_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes_open, name = "Sex") +
scale_x_discrete(labels = my_labels_full) +
ylim(19, 61) +
savs_ggplot_theme -> full_mass_emmean_plot
full_mass_emmean_plot
Mass - F
ggplot() +
geom_jitter(data = dat_F,
aes(x = month,
y = mass_g,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_F[full_emmeans_F$measurement == "Mass",],
aes(x = month,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_F[full_emmeans_F$measurement == "Mass",],
aes(x = month,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab("Lizard Mass (g)") +
# sex differences
annotate(geom = "text", x = 3, y = 54, label = "*", size = 6) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 61, label = "a", size = 3) +
annotate(geom = "text", x = 2, y = 61, label = "b", size = 3) +
annotate(geom = "text", x = 3, y = 61, label = "a", size = 3) +
scale_color_manual(values = Fem, name = "Sex") +
scale_fill_manual(values = Fem, name = "Sex") +
scale_shape_manual(values = 19, name = "Sex") +
scale_x_discrete(labels = my_labels_full_F) +
ylim(19, 61) +
savs_ggplot_theme -> full_mass_emmean_plot_F
full_mass_emmean_plot_F
Mass - M
ggplot() +
geom_jitter(data = dat_M,
aes(x = month,
y = mass_g,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_M[full_emmeans_M$measurement == "Mass",],
aes(x = month,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_M[full_emmeans_M$measurement == "Mass",],
aes(x = month,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab("Lizard Mass (g)") +
# sex differences
annotate(geom = "text", x = 3, y = 54, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 1, y = 61, label = "a,b", size = 3) +
annotate(geom = "text", x = 2, y = 61, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 61, label = "b", size = 3) +
scale_color_manual(values = Male, name = "Sex") +
scale_fill_manual(values = Male, name = "Sex") +
scale_shape_manual(values = 15, name = "Sex") +
scale_x_discrete(labels = my_labels_full_M) +
ylim(19, 61) +
savs_ggplot_theme -> full_mass_emmean_plot_M
full_mass_emmean_plot_M
Body Condition
ggplot() +
geom_jitter(data = dat,
aes(x = month_sex,
y = SMI,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_all[full_emmeans_all$measurement == "SMI",],
aes(x = month_sex,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_all[full_emmeans_all$measurement == "SMI",],
aes(x = month_sex,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab("Body Condition (g')") +
# sex differences
annotate(geom = "text", x = 1.5, y = 44, label = "*", size = 6) +
annotate(geom = "text", x = 3.5, y = 44, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 2, y = 49, label = "A", size = 3) +
annotate(geom = "text", x = 4, y = 49, label = "A", size = 3) +
annotate(geom = "text", x = 6, y = 49, label = "B", size = 3) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 49, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 49, label = "b", size = 3) +
annotate(geom = "text", x = 5, y = 49, label = "a", size = 3) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_fill_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes_open, name = "Sex") +
scale_x_discrete(labels = my_labels_full) +
ylim(19, 50) +
savs_ggplot_theme -> full_SMI_emmean_plot
full_SMI_emmean_plot
Body Cond - F
ggplot() +
geom_jitter(data = dat_F,
aes(x = month,
y = SMI,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_F[full_emmeans_F$measurement == "SMI",],
aes(x = month,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_F[full_emmeans_F$measurement == "SMI",],
aes(x = month,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab("Body Condition (g')") +
# sex differences
annotate(geom = "text", x = 1, y = 47, label = "*", size = 6) +
annotate(geom = "text", x = 2, y = 47, label = "*", size = 6) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 50, label = "a", size = 3) +
annotate(geom = "text", x = 2, y = 50, label = "b", size = 3) +
annotate(geom = "text", x = 3, y = 50, label = "a", size = 3) +
scale_color_manual(values = Fem, name = NULL) +
scale_fill_manual(values = Fem, name = NULL) +
scale_shape_manual(values = 21, name = NULL) +
scale_x_discrete(labels = my_labels_full_F) +
ylim(19, 50) +
savs_ggplot_theme +
theme(plot.margin = margin(t = 0, r = 0, b = 0, l = 10, unit = "pt"),
legend.margin = margin(20,0,0,0, unit = "pt")) -> full_SMI_emmean_plot_F
full_SMI_emmean_plot_F
Body Cond - M
ggplot() +
geom_jitter(data = dat_M,
aes(x = month,
y = SMI,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_M[full_emmeans_M$measurement == "SMI",],
aes(x = month,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_M[full_emmeans_M$measurement == "SMI",],
aes(x = month,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab("Body Condition (g')") +
# sex differences
annotate(geom = "text", x = 1, y = 47, label = "*", size = 6) +
annotate(geom = "text", x = 2, y = 47, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 1, y = 50, label = "a", size = 3) +
annotate(geom = "text", x = 2, y = 50, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 50, label = "b", size = 3) +
scale_color_manual(values = Male, name = NULL) +
scale_fill_manual(values = Male, name = NULL) +
scale_shape_manual(values = 22, name = NULL) +
scale_x_discrete(labels = my_labels_full_M) +
ylim(19, 50) +
savs_ggplot_theme +
theme(plot.margin = margin(t = 0, r = 0, b = 0, l = 0, unit = "pt"),
legend.margin = margin(20,0,0,0, unit = "pt")) -> full_SMI_emmean_plot_M
full_SMI_emmean_plot_M
Osmolality
ggplot() +
geom_jitter(data = dat,
aes(x = month_sex,
y = osmolality_mmol_kg_mean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_all[full_emmeans_all$measurement == "Osmolality",],
aes(x = month_sex,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_all[full_emmeans_all$measurement == "Osmolality",],
aes(x = month_sex,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab(bquote('Osmolality (mmol '*kg^-1*')')) +
# sex differences
annotate(geom = "text", x = 1.5, y = 422, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 2, y = 445, label = "AB", size = 3) +
annotate(geom = "text", x = 4, y = 445, label = "A", size = 3) +
annotate(geom = "text", x = 6, y = 445, label = "B", size = 3) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 445, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 445, label = "a", size = 3) +
annotate(geom = "text", x = 5, y = 445, label = "a", size = 3) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_fill_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes_open, name = "Sex") +
scale_x_discrete(labels = my_labels_full) +
scale_y_continuous(limits = c(320, 450), breaks = c(320, 360, 400, 440)) +
savs_ggplot_theme -> full_osml_emmean_plot
full_osml_emmean_plot## Warning: Removed 5 rows containing missing values (`geom_point()`).
Osmolality - F
ggplot() +
geom_jitter(data = dat_F,
aes(x = month,
y = osmolality_mmol_kg_mean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_F[full_emmeans_F$measurement == "Osmolality",],
aes(x = month,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_F[full_emmeans_F$measurement == "Osmolality",],
aes(x = month,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab(bquote('Osmolality (mmol '*kg^-1*')')) +
# sex differences
annotate(geom = "text", x = 1, y = 435, label = "*", size = 6) +
# time differences for FEMALES
annotate(geom = "text", x = 1, y = 448, label = "a", size = 3) +
annotate(geom = "text", x = 2, y = 448, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 448, label = "a", size = 3) +
scale_color_manual(values = Fem, name = "Sex") +
scale_fill_manual(values = Fem, name = "Sex") +
scale_shape_manual(values = 21, name = "Sex") +
scale_x_discrete(labels = my_labels_full_F) +
scale_y_continuous(limits = c(320, 450), breaks = c(320, 360, 400, 440)) +
savs_ggplot_theme +
theme(plot.margin = margin(t = 0, r = 0, b = 0, l = 0, unit = "pt"),
legend.margin = margin(0,0,0,0, unit = "pt")) -> full_osml_emmean_plot_F
full_osml_emmean_plot_F## Warning: Removed 3 rows containing missing values (`geom_point()`).
ggsave(filename = "szn_change_osml_F_standalone.pdf",
plot = full_osml_emmean_plot_F,
path = "./results_figures",
device = "pdf",
dpi = 600,
units = "mm",
width = 60, height = 80
)## Warning: Removed 3 rows containing missing values (`geom_point()`).Osmolality - M
ggplot() +
geom_jitter(data = dat_M,
aes(x = month,
y = osmolality_mmol_kg_mean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_M[full_emmeans_M$measurement == "Osmolality",],
aes(x = month,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_M[full_emmeans_M$measurement == "Osmolality",],
aes(x = month,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
xlab("") +
ylab(bquote('Osmolality (mmol '*kg^-1*')')) +
# sex differences
annotate(geom = "text", x = 1, y = 435, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 1, y = 448, label = "a,b", size = 3) +
annotate(geom = "text", x = 2, y = 448, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 448, label = "b", size = 3) +
scale_color_manual(values = Male, name = "Sex") +
scale_fill_manual(values = Male, name = "Sex") +
scale_shape_manual(values = 22, name = "Sex") +
scale_x_discrete(labels = my_labels_full_M) +
scale_y_continuous(limits = c(320, 450), breaks = c(320, 360, 400, 440)) +
savs_ggplot_theme +
theme(plot.margin = margin(t = 0, r = 0, b = 0, l = 0, unit = "pt"),
legend.margin = margin(0,0,0,0, unit = "pt")) -> full_osml_emmean_plot_M
full_osml_emmean_plot_M## Warning: Removed 2 rows containing missing values (`geom_point()`).
ggsave(filename = "szn_change_osml_M_standalone.pdf",
plot = full_osml_emmean_plot_M,
path = "./results_figures",
device = "pdf",
dpi = 600,
units = "mm",
width = 60, height = 80
)## Warning: Removed 2 rows containing missing values (`geom_point()`).CEWL
ggplot() +
geom_jitter(data = dat,
aes(x = month_sex,
y = CEWL_g_m2h,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_all[full_emmeans_all$measurement == "CEWL",],
aes(x = month_sex,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_all[full_emmeans_all$measurement == "CEWL",],
aes(x = month_sex,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
# sex differences
annotate(geom = "text", x = 3.5, y = 20, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 2, y = 24, label = "A", size = 3) +
annotate(geom = "text", x = 4, y = 24, label = "B", size = 3) +
annotate(geom = "text", x = 6, y = 24, label = "A", size = 3) +
# time differences for FEMALES = none
annotate(geom = "text", x = 1, y = 24, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 24, label = "a", size = 3) +
annotate(geom = "text", x = 5, y = 24, label = "a", size = 3) +
xlab("") +
ylab(bquote('CEWL (g '*m^-2*' '*h^-1*')')) +
ylim(0, 25) +
scale_color_manual(values = my_colors, name = "Sex") +
scale_fill_manual(values = my_colors, name = "Sex") +
scale_shape_manual(values = my_shapes_open, name = "Sex") +
scale_x_discrete(labels = my_labels_full) +
savs_ggplot_theme -> full_CEWL_emmean_plot
full_CEWL_emmean_plot## Warning: Removed 1 rows containing missing values (`geom_point()`).
CEWL - F
ggplot() +
geom_jitter(data = dat_F,
aes(x = month,
y = CEWL_g_m2h,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_F[full_emmeans_F$measurement == "CEWL",],
aes(x = month,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_F[full_emmeans_F$measurement == "CEWL",],
aes(x = month,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
# sex differences
annotate(geom = "text", x = 2, y = 22.5, label = "*", size = 6) +
# time differences for FEMALES = none
annotate(geom = "text", x = 1, y = 25, label = "a", size = 3) +
annotate(geom = "text", x = 2, y = 25, label = "a", size = 3) +
annotate(geom = "text", x = 3, y = 25, label = "a", size = 3) +
xlab("") +
ylab(bquote('CEWL (g '*m^-2*' '*h^-1*')')) +
ylim(0, 25) +
scale_color_manual(values = Fem, name = "Sex") +
scale_fill_manual(values = Fem, name = "Sex") +
scale_shape_manual(values = 21, name = "Sex") +
scale_x_discrete(labels = my_labels_full_F) +
savs_ggplot_theme +
theme(plot.margin = margin(t = 0, r = 0, b = 0, l = 4.5, unit = "pt"),
legend.margin = margin(0,0,0,0, unit = "pt")) -> full_CEWL_emmean_plot_F
full_CEWL_emmean_plot_F
CEWL - M
ggplot() +
geom_jitter(data = dat_M,
aes(x = month,
y = CEWL_g_m2h,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 0.3,
size = 0.9,
position = position_jitter(height = 0, width = 0.2)) +
geom_point(data = full_emmeans_M[full_emmeans_M$measurement == "CEWL",],
aes(x = month,
y = emmean,
shape = sex_M_F,
color = sex_M_F,
fill = sex_M_F),
alpha = 1,
size = 2) +
geom_errorbar(data = full_emmeans_M[full_emmeans_M$measurement == "CEWL",],
aes(x = month,
y = emmean,
ymin = lower.CL,
ymax = upper.CL,
color = sex_M_F),
width = .4) +
# sex differences
annotate(geom = "text", x = 2, y = 22.5, label = "*", size = 6) +
# time differences for MALES
annotate(geom = "text", x = 1, y = 25, label = "a", size = 3) +
annotate(geom = "text", x = 2, y = 25, label = "b", size = 3) +
annotate(geom = "text", x = 3, y = 25, label = "a", size = 3) +
xlab("") +
ylab(bquote('CEWL (g '*m^-2*' '*h^-1*')')) +
ylim(0, 25) +
scale_color_manual(values = Male, name = "Sex") +
scale_fill_manual(values = Male, name = "Sex") +
scale_shape_manual(values = 22, name = "Sex") +
scale_x_discrete(labels = my_labels_full_M) +
savs_ggplot_theme +
theme(plot.margin = margin(t = 0, r = 0, b = 0, l = 0, unit = "pt"),
legend.margin = margin(0,0,0,0, unit = "pt")) -> full_CEWL_emmean_plot_M
full_CEWL_emmean_plot_M## Warning: Removed 1 rows containing missing values (`geom_point()`).
ggsave(filename = "szn_change_CEWL_M_standalone.pdf",
plot = full_CEWL_emmean_plot_M,
path = "./results_figures",
device = "pdf",
dpi = 600,
units = "mm",
width = 60, height = 80
)## Warning: Removed 1 rows containing missing values (`geom_point()`).Arrange and Export M&F Separate Panel Plot
# CEWL
full_CEWL_emmean_plot_M_multi <- full_CEWL_emmean_plot_M +
theme(axis.text.y = element_blank(),
axis.ticks.y = element_blank(),
axis.text.x = element_blank(),
axis.ticks.x = element_blank(),
legend.position = "none") +
labs(y = NULL)
full_CEWL_emmean_plot_F_multi <- full_CEWL_emmean_plot_F +
theme(axis.text.x = element_blank(),
axis.ticks.x = element_blank(),
legend.position = "none")
# osmolality
full_osml_emmean_plot_M_multi <- full_osml_emmean_plot_M +
theme(axis.text.y = element_blank(),
axis.ticks.y = element_blank(),
axis.text.x = element_blank(),
axis.ticks.x = element_blank(),
legend.position = "none") +
labs(y = NULL)
full_osml_emmean_plot_F_multi <- full_osml_emmean_plot_F +
theme(axis.text.x = element_blank(),
axis.ticks.x = element_blank(),
legend.position = "none")
# body condition
full_SMI_emmean_plot_M_multi <- full_SMI_emmean_plot_M +
theme(axis.text.y = element_blank(),
axis.ticks.y = element_blank()) +
labs(y = NULL)
# arrange M & F together
ggarrange(full_CEWL_emmean_plot_F_multi, full_CEWL_emmean_plot_M_multi,
full_osml_emmean_plot_F_multi, full_osml_emmean_plot_M_multi,
full_SMI_emmean_plot_F, full_SMI_emmean_plot_M_multi,
ncol = 2, nrow = 3,
widths = c(1.36, 1), # make x-axis widths match when exported
heights = c(1, 1, 1.36), # heights
labels = c("(a)", "", "(b)", "", "(c)", ""),
label.y = 1.02,
label.x = -0.05,
font.label = list(color = "black", face = "bold", size = 10)
) -> full_szn_change_multi_fig_best## Warning: Removed 1 rows containing missing values (`geom_point()`).## Warning: Removed 3 rows containing missing values (`geom_point()`).## Warning: Removed 2 rows containing missing values (`geom_point()`).Export Full Data Multi-Panel Fig
ggarrange(full_CEWL_emmean_plot,
full_osml_emmean_plot,
#full_mass_emmean_plot,
full_SMI_emmean_plot,
ncol = 3, nrow = 1,
#ncol = 1, nrow = 3,
labels = c("A", "B", "C"),
common.legend = TRUE,
legend = "bottom"
) -> full_szn_change_multi_fig_h## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Removed 1 rows containing missing values (`geom_point()`).## Warning: Removed 5 rows containing missing values (`geom_point()`).Export Full & Repeats Multi-Panel Fig
ggarrange(full_CEWL_emmean_plot,
full_osml_emmean_plot,
full_mass_emmean_plot,
full_SMI_emmean_plot,
CEWL_emmean_plot,
osml_emmean_plot,
mass_emmean_plot,
SMI_emmean_plot,
ncol = 4, nrow = 2,
# ncol = 1, nrow = 3,
labels = c("full", "full", "full", "full",
"rep", "rep", "rep", "rep"),
common.legend = TRUE,
legend = "bottom"
) -> full_vs_rep_szn_change_multi_fig## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Removed 1 rows containing missing values (`geom_point()`).## Warning: Removed 5 rows containing missing values (`geom_point()`).## Warning: Removed 1 rows containing missing values (`geom_point()`).## Warning: Removed 2 rows containing missing values (`geom_point()`).Reproductive Effort
Is reproductive effort related to water balance?
Egg Investment
## capture_date_time capture_date radio_collar_serial
## Min. :2021-04-23 04:39:00.00 Min. :2021-04-23 Length:18
## 1st Qu.:2021-04-23 08:07:30.00 1st Qu.:2021-04-23 Class :character
## Median :2021-04-24 05:51:30.00 Median :2021-05-01 Mode :character
## Mean :2021-04-29 23:59:56.25 Mean :2021-04-30
## 3rd Qu.:2021-05-08 02:32:30.00 3rd Qu.:2021-05-08
## Max. :2021-05-08 05:27:00.00 Max. :2021-05-08
## NA's :2
## PIT_tag_ID individual_ID sex_M_F mass_g SVL_mm
## Length:18 F-01 :2 Female:18 Min. :31.00 Min. : 90.0
## Class :character F-03 :2 Male : 0 1st Qu.:35.60 1st Qu.: 95.0
## Mode :character F-05 :2 Median :40.00 Median :100.0
## F-06 :2 Mean :41.73 Mean :100.7
## F-11 :2 3rd Qu.:46.75 3rd Qu.:106.5
## F-12 :2 Max. :58.80 Max. :114.0
## (Other):6
## tmt tmt_mass_change_g capture_to_msmt hematocrit_percent
## Water Tmt: 8 Min. :-3.0000 Min. : 25.00 Min. :23.00
## Sham Tmt :10 1st Qu.:-1.4000 1st Qu.: 55.25 1st Qu.:30.25
## No Tmt : 0 Median :-0.1000 Median : 147.00 Median :32.50
## Mean :-0.5941 Mean : 225.88 Mean :33.78
## 3rd Qu.: 0.0000 3rd Qu.: 247.00 3rd Qu.:37.00
## Max. : 1.8000 Max. :1638.00 Max. :44.00
## NA's :1 NA's :2
## cloacal_temp_C CEWL_g_m2h msmt_temp_C msmt_RH_percent
## Min. :28.00 Min. : 1.493 Min. :23.13 Min. :12.20
## 1st Qu.:33.25 1st Qu.: 8.986 1st Qu.:30.47 1st Qu.:14.19
## Median :34.00 Median :10.320 Median :31.59 Median :14.90
## Mean :33.72 Mean :11.109 Mean :30.79 Mean :16.29
## 3rd Qu.:35.00 3rd Qu.:13.559 3rd Qu.:31.89 3rd Qu.:18.02
## Max. :36.00 Max. :21.673 Max. :32.23 Max. :28.17
##
## e_s_kPa msmt_VPD_kPa osmolality_mmol_kg_mean month week
## Min. :2.830 Min. :2.033 Min. :329.5 April:9 1 :9
## 1st Qu.:4.358 1st Qu.:3.692 1st Qu.:338.0 May :9 3 :9
## Median :4.646 Median :3.916 Median :345.0 July :0 12:0
## Mean :4.463 Mean :3.751 Mean :350.9
## 3rd Qu.:4.725 3rd Qu.:4.054 3rd Qu.:362.9
## Max. :4.818 Max. :4.167 Max. :390.0
##
## week_num ultrasound_date gravid_Y_N number_eggs
## Min. :1 Min. :2021-04-24 Not Gravid: 6 Min. :3.000
## 1st Qu.:1 1st Qu.:2021-04-24 Gravid :12 1st Qu.:3.750
## Median :2 Median :2021-05-01 Median :4.000
## Mean :2 Mean :2021-05-01 Mean :3.833
## 3rd Qu.:3 3rd Qu.:2021-05-08 3rd Qu.:4.000
## Max. :3 Max. :2021-05-08 Max. :5.000
## NA's :6
## largest_egg_diameter_cm largest_egg_circumference_cm stage
## Min. :0.640 Min. :1.820 small round:1
## 1st Qu.:0.865 1st Qu.:2.542 large round:6
## Median :1.095 Median :3.185 soft :2
## Mean :1.173 Mean :3.292 soft oblong:1
## 3rd Qu.:1.488 3rd Qu.:4.040 firm oblong:1
## Max. :1.800 Max. :4.740 hard oblong:1
## NA's :6 NA's :6 NA's :6
## dev_point SMI month_sex date_chunk_pretty
## Min. :1.00 Min. :29.24 April Female:9 Apr 23-25:9
## 1st Qu.:2.00 1st Qu.:34.63 April Male :0 May 7-9 :9
## Median :2.00 Median :36.73 May Female :9 Jul 14 :0
## Mean :2.75 Mean :37.83 May Male :0
## 3rd Qu.:3.25 3rd Qu.:41.85 July Female :0
## Max. :5.00 Max. :45.63 July Male :0
## NA's :6
## percent_water_mass_change n
## Min. :-6.834 2:18
## 1st Qu.:-4.242 1: 0
## Median :-0.250
## Mean :-1.539
## 3rd Qu.: 0.000
## Max. : 3.462
## NA's :1dat_repeat_females_eggs <- dat_repeat_females %>%
dplyr::filter(gravid_Y_N == "Gravid") %>%
dplyr::filter(month == "May") # so we don't have any repeats
summary(lm(data = dat_repeat_females_eggs, number_eggs ~ mass_g))##
## Call:
## lm(formula = number_eggs ~ mass_g, data = dat_repeat_females_eggs)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.2855 -0.2647 -0.1891 0.1467 0.7541
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.11895 0.69935 1.600 0.1607
## mass_g 0.05975 0.01554 3.846 0.0085 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4103 on 6 degrees of freedom
## Multiple R-squared: 0.7114, Adjusted R-squared: 0.6633
## F-statistic: 14.79 on 1 and 6 DF, p-value: 0.008501##
## Call:
## lm(formula = number_eggs ~ SVL_mm, data = dat_repeat_females_eggs)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.38064 -0.30282 -0.00921 0.28214 0.44492
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -5.43434 2.07306 -2.621 0.03951 *
## SVL_mm 0.09173 0.02066 4.439 0.00438 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.369 on 6 degrees of freedom
## Multiple R-squared: 0.7666, Adjusted R-squared: 0.7277
## F-statistic: 19.71 on 1 and 6 DF, p-value: 0.00438##
## Call:
## lm(formula = number_eggs ~ SMI, data = dat_repeat_females_eggs)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.9073 -0.3127 -0.0824 0.2596 0.9036
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.04386 2.40577 0.018 0.986
## SMI 0.09238 0.05970 1.548 0.173
##
## Residual standard error: 0.6457 on 6 degrees of freedom
## Multiple R-squared: 0.2853, Adjusted R-squared: 0.1661
## F-statistic: 2.395 on 1 and 6 DF, p-value: 0.1727## Analysis of Variance Table
##
## Response: number_eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## mass_g 1 2.4899 2.48991 14.79 0.008501 **
## Residuals 6 1.0101 0.16835
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## Analysis of Variance Table
##
## Response: number_eggs
## Df Sum Sq Mean Sq F value Pr(>F)
## SVL_mm 1 2.68307 2.68307 19.706 0.00438 **
## Residuals 6 0.81693 0.13616
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Repeat Individuals
I have 9 females measured at both April and May time points. Let’s see if their physiology changes consistently based on becoming gravid.
First, calculate the change in CEWL, osmolality, and mass for each of those females between April to May.
fem_starting_values <- dat_repeat_females %>%
# get the values for April measurements only
dplyr::filter(month == "April") %>%
dplyr::select(individual_ID,
CEWL_April = CEWL_g_m2h,
osml_April = osmolality_mmol_kg_mean,
mass_April = mass_g,
SMI_April = SMI)
fem_change <- dat_repeat_females %>%
# add starting values
left_join(fem_starting_values, by = 'individual_ID') %>%
# make sure we have complete data
dplyr::filter(complete.cases(mass_April, CEWL_April)) %>%
# calculate difference for each msmt versus in April
mutate(mass_change_from_April = mass_g - mass_April,
SMI_change_from_April = SMI - SMI_April,
CEWL_change_from_April = CEWL_g_m2h - CEWL_April,
osml_change_from_April = osmolality_mmol_kg_mean - osml_April)
summary(fem_change)## capture_date_time capture_date radio_collar_serial
## Min. :2021-04-23 04:39:00.00 Min. :2021-04-23 Length:18
## 1st Qu.:2021-04-23 08:07:30.00 1st Qu.:2021-04-23 Class :character
## Median :2021-04-24 05:51:30.00 Median :2021-05-01 Mode :character
## Mean :2021-04-29 23:59:56.25 Mean :2021-04-30
## 3rd Qu.:2021-05-08 02:32:30.00 3rd Qu.:2021-05-08
## Max. :2021-05-08 05:27:00.00 Max. :2021-05-08
## NA's :2
## PIT_tag_ID individual_ID sex_M_F mass_g SVL_mm
## Length:18 F-01 :2 Female:18 Min. :31.00 Min. : 90.0
## Class :character F-03 :2 Male : 0 1st Qu.:35.60 1st Qu.: 95.0
## Mode :character F-05 :2 Median :40.00 Median :100.0
## F-06 :2 Mean :41.73 Mean :100.7
## F-11 :2 3rd Qu.:46.75 3rd Qu.:106.5
## F-12 :2 Max. :58.80 Max. :114.0
## (Other):6
## tmt tmt_mass_change_g capture_to_msmt hematocrit_percent
## Water Tmt: 8 Min. :-3.0000 Min. : 25.00 Min. :23.00
## Sham Tmt :10 1st Qu.:-1.4000 1st Qu.: 55.25 1st Qu.:30.25
## No Tmt : 0 Median :-0.1000 Median : 147.00 Median :32.50
## Mean :-0.5941 Mean : 225.88 Mean :33.78
## 3rd Qu.: 0.0000 3rd Qu.: 247.00 3rd Qu.:37.00
## Max. : 1.8000 Max. :1638.00 Max. :44.00
## NA's :1 NA's :2
## cloacal_temp_C CEWL_g_m2h msmt_temp_C msmt_RH_percent
## Min. :28.00 Min. : 1.493 Min. :23.13 Min. :12.20
## 1st Qu.:33.25 1st Qu.: 8.986 1st Qu.:30.47 1st Qu.:14.19
## Median :34.00 Median :10.320 Median :31.59 Median :14.90
## Mean :33.72 Mean :11.109 Mean :30.79 Mean :16.29
## 3rd Qu.:35.00 3rd Qu.:13.559 3rd Qu.:31.89 3rd Qu.:18.02
## Max. :36.00 Max. :21.673 Max. :32.23 Max. :28.17
##
## e_s_kPa msmt_VPD_kPa osmolality_mmol_kg_mean month week
## Min. :2.830 Min. :2.033 Min. :329.5 April:9 1 :9
## 1st Qu.:4.358 1st Qu.:3.692 1st Qu.:338.0 May :9 3 :9
## Median :4.646 Median :3.916 Median :345.0 July :0 12:0
## Mean :4.463 Mean :3.751 Mean :350.9
## 3rd Qu.:4.725 3rd Qu.:4.054 3rd Qu.:362.9
## Max. :4.818 Max. :4.167 Max. :390.0
##
## week_num ultrasound_date gravid_Y_N number_eggs
## Min. :1 Min. :2021-04-24 Not Gravid: 6 Min. :3.000
## 1st Qu.:1 1st Qu.:2021-04-24 Gravid :12 1st Qu.:3.750
## Median :2 Median :2021-05-01 Median :4.000
## Mean :2 Mean :2021-05-01 Mean :3.833
## 3rd Qu.:3 3rd Qu.:2021-05-08 3rd Qu.:4.000
## Max. :3 Max. :2021-05-08 Max. :5.000
## NA's :6
## largest_egg_diameter_cm largest_egg_circumference_cm stage
## Min. :0.640 Min. :1.820 small round:1
## 1st Qu.:0.865 1st Qu.:2.542 large round:6
## Median :1.095 Median :3.185 soft :2
## Mean :1.173 Mean :3.292 soft oblong:1
## 3rd Qu.:1.488 3rd Qu.:4.040 firm oblong:1
## Max. :1.800 Max. :4.740 hard oblong:1
## NA's :6 NA's :6 NA's :6
## dev_point SMI month_sex date_chunk_pretty
## Min. :1.00 Min. :29.24 April Female:9 Apr 23-25:9
## 1st Qu.:2.00 1st Qu.:34.63 April Male :0 May 7-9 :9
## Median :2.00 Median :36.73 May Female :9 Jul 14 :0
## Mean :2.75 Mean :37.83 May Male :0
## 3rd Qu.:3.25 3rd Qu.:41.85 July Female :0
## Max. :5.00 Max. :45.63 July Male :0
## NA's :6
## percent_water_mass_change n CEWL_April osml_April
## Min. :-6.834 2:18 Min. : 1.493 Min. :329.5
## 1st Qu.:-4.242 1: 0 1st Qu.: 8.940 1st Qu.:333.3
## Median :-0.250 Median :10.353 Median :338.0
## Mean :-1.539 Mean :11.304 Mean :344.8
## 3rd Qu.: 0.000 3rd Qu.:13.640 3rd Qu.:347.7
## Max. : 3.462 Max. :21.673 Max. :390.0
## NA's :1
## mass_April SMI_April mass_change_from_April SMI_change_from_April
## Min. :31.00 Min. :29.24 Min. :-2.100 Min. :-0.6421
## 1st Qu.:33.00 1st Qu.:32.19 1st Qu.: 0.000 1st Qu.: 0.0000
## Median :40.00 Median :35.13 Median : 0.000 Median : 0.0000
## Mean :39.44 Mean :34.93 Mean : 2.289 Mean : 2.8985
## 3rd Qu.:44.00 3rd Qu.:37.02 3rd Qu.: 4.425 3rd Qu.: 3.8356
## Max. :52.00 Max. :41.30 Max. : 9.000 Max. :14.2653
##
## CEWL_change_from_April osml_change_from_April
## Min. :-8.2493 Min. :-50.667
## 1st Qu.:-1.5265 1st Qu.: 0.000
## Median : 0.0000 Median : 0.000
## Mean :-0.1947 Mean : 6.111
## 3rd Qu.: 0.0000 3rd Qu.: 11.625
## Max. : 7.8247 Max. : 45.000
## Plot Indiv Change
ggplot(fem_change) +
aes(x = month,
y = mass_change_from_April,
group = individual_ID,
color = gravid_Y_N
) +
geom_point() +
geom_line() +
theme_bw()
ggplot(fem_change) +
aes(x = month,
y = CEWL_change_from_April,
group = individual_ID,
color = gravid_Y_N
) +
geom_point() +
geom_line() +
theme_bw()
ggplot(fem_change) +
aes(x = month,
y = osml_change_from_April,
group = individual_ID,
color = gravid_Y_N
) +
geom_point() +
geom_line() +
theme_bw()
On average, mass and osmolality both increase. There’s only one lizard who was not gravid at both time points, so we can’t compare change based on becoming gravid or not. I think it is reasonable to conclude that on average, gravid females gained mass and became less hydrated. But, we cannot say that this was DUE to gravidity, because we have nothing to compare it to
Stats Indiv Change
Run stats tests on the change rates for only the females who were gravid in May, and use only the actual delta values, in their rows for May data.
# check osml outlier range
fem_change %>%
group_by(month) %>%
summarise(mean(osml_change_from_April), sd(osml_change_from_April))## # A tibble: 2 × 3
## month `mean(osml_change_from_April)` `sd(osml_change_from_April)`
## <fct> <dbl> <dbl>
## 1 April 0 0
## 2 May 12.2 29.0# for stats
fem_change_stats <- fem_change %>%
# don't use the one female still non-gravid in May
dplyr::filter(individual_ID != "F-12") %>%
# only use the >resulting< change values
dplyr::filter(month == "May") %>%
# remove the crazy low outlier based on sd calculated above
dplyr::filter(osml_change_from_April > (12-(29*2)))ttests:
##
## One Sample t-test
##
## data: fem_change_stats$CEWL_change_from_April
## t = -0.59229, df = 6, p-value = 0.5753
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## -5.873479 3.584194
## sample estimates:
## mean of x
## -1.144643##
## One Sample t-test
##
## data: fem_change_stats$osml_change_from_April
## t = 2.8997, df = 6, p-value = 0.02734
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## 3.28688 38.80836
## sample estimates:
## mean of x
## 21.04762##
## One Sample t-test
##
## data: fem_change_stats$mass_change_from_April
## t = 3.5281, df = 6, p-value = 0.0124
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## 1.514744 8.370970
## sample estimates:
## mean of x
## 4.942857##
## One Sample t-test
##
## data: fem_change_stats$SMI_change_from_April
## t = 2.7067, df = 6, p-value = 0.03526
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## 0.5695776 11.2980848
## sample estimates:
## mean of x
## 5.933831Stats Likelihood
Does anything correlate with the likelihood of being/becoming gravid?
Start with trying each possible logistic model to predict gravidity:
repro_logit1 <- lme4::glmer(data = dat_repro_females,
gravid_Y_N ~ mass_g +
(1|individual_ID),
family = "binomial")
repro_logit2 <- lme4::glmer(data = dat_repro_females,
gravid_Y_N ~ SMI +
(1|individual_ID),
family = "binomial")
repro_logit3 <- lme4::glmer(data = dat_repro_females,
gravid_Y_N ~ osmolality_mmol_kg_mean +
(1|individual_ID),
family = "binomial")## boundary (singular) fit: see help('isSingular')repro_logit4 <- lme4::glmer(data = dat_repro_females,
gravid_Y_N ~ CEWL_g_m2h +
(1|individual_ID),
family = "binomial")
car::Anova(repro_logit1, type = 2,
test.statistic = "Chisq", error.estimate = "pearson")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: gravid_Y_N
## Chisq Df Pr(>Chisq)
## mass_g 2.3158 1 0.1281## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: gravid_Y_N
## Chisq Df Pr(>Chisq)
## SMI 1.3948 1 0.2376## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: gravid_Y_N
## Chisq Df Pr(>Chisq)
## osmolality_mmol_kg_mean 0.2219 1 0.6376## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: gravid_Y_N
## Chisq Df Pr(>Chisq)
## CEWL_g_m2h 0.0045 1 0.9468What about an ANOVA?
Plot Likelihood Gravid ~ Mass
dat_repro_females %>%
# dummy-fy gravid variable
mutate(gravid_Y_N = factor(gravid_Y_N, levels = c("Not Gravid", "Gravid"),
labels = c(0,1)),
# goes 1-2, so fix to be 0-1
gravid_Y_N = as.numeric(gravid_Y_N)-1) %>%
# nooooowwwww I can plot and the stat_smooth will actually work
ggplot() +
aes(x = mass_g,
y = gravid_Y_N) +
geom_jitter(alpha=.5,
position = position_jitter(height = 0.1, width = 0)
) +
# must make response variable 0-1 in order for this to work
stat_smooth(geom = "smooth",
method = "glm",
formula = y~x,
se = FALSE,
method.args = list(family = "quasibinomial")
) +
ylab("") +
scale_y_continuous(labels = c("Not Gravid", "", "", "", "Gravid"),
limits = c(0,1)) +
xlab("Body Mass (g)") +
savs_ggplot_theme## Warning: Removed 11 rows containing missing values (`geom_point()`).
Summary Figure
dodger <- position_dodge(width = 0.2) # make dodge amt for both points and lines
dat_repro_females %>%
# only use data for April-May
dplyr::filter(month %in% c("April", "May")) %>%
ggplot() +
aes(x = month,
y = mass_g,
color = gravid_Y_N,
shape = gravid_Y_N,
alpha = gravid_Y_N,
group = individual_ID) +
geom_line(position = dodger) +
geom_point(size = 2,
position = dodger, #position_jitter(height = 0.01, width = 0),
# alpha = 0.5
) +
labs(y = "Body Mass (g)",
x = NULL) +
scale_color_manual(values = c(Fem_nongrav, Fem_grav), name = NULL) +
scale_shape_manual(values = c(19, 19), name = NULL) +
scale_alpha_manual(values = c(0.3, 0.7), name = NULL) +
savs_ggplot_theme +
theme(legend.position = "right",
plot.margin = margin(t = 1.5, r = 0.5, b = 0.5, l = 0.5, unit = "pt"),
legend.margin = margin(0,0,0,0, unit = "pt")) -> when_gravid
when_gravid